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Figure 1. USAF Firefighters Drill. 






The autonomous firefighter manipulation problem has two main areas: design 
of a manipulator mechanism and controller design for that mechanism. In 
addressing these problem areas the following assumptions were made: 

1 . The manipulator must manipulate the USAF tool shown in Figure 1 . 

2. The manipulator must be as autonomous as possible. 

3. Due to (2) , the maniplator must be as power efficient as possible in 
order to conserve on-board power supplies. 

4. The tool positioning requirement would be approximately + 0.10 in and 
+ 5 deg at the drill tip. 

The mechanism design must first be concerned with the expected loads 
associated with the drilling process. Since these are, in general, widely 
variable depending on the type of drill, type of material to be drilled, 
incidence angle, reinforcements, and many other parameters, we were not 
surprized to find very litle published material on the drilling process in the 
open literature. To remedy this we undertook several experimental projects to 
measure and characterize the particular process of drilling aircraft skins. 
Given the measured drilling load characteristics, the manipulator mechanism 
and actuators could be designed or selected. 

The manipulator controller is concerned with guiding the manipulator from 
point-to-point in a power efficient manner, and maintaining a tool orientation 
while a hole is drilled. The following discussion is about the first part of 
the controller problem, the point-to-point control and the tasks associated 
with solving this problem. Clearly, the manipulation requirements outlined in 
assumption 4 above could not be considered high performance specifications. 
Instead, the controller of interest was one aimed at effective coordination of 
movement and was therefore called a low performance controller. 
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The remainder of this report is divided into three sections: the next 
section summarizes the major results acheived to data; the following section 
describes the most recent mainpulator controller development project; the last 
section discusses final drilling load measurement studies. 

II. Results 

A. Facilities Developed 

Perhaps the most important facility that has been developed for this 
project is the Low Performance Control (LPC) Lab. This test bed facility 
consists of a computer based controller system designed to be used to 
investigate the merits of various candidate controller and sensor 
philosophies for low performance plants. The next section of this report 
describes the LPC lab in more detail when coupled to a manipulator so it 
will only be mentioned briefly here. The Lab presently consists of an IBM 
PC- AT connected to a rigid, revolute robot (Neptune by Alan Hayes Corp.). 
The Computer /manipulator interaction is accomplished via LABPAC hardware 
and software. In addition, hardware for optical encoding of joint data, 
and an improved data taking and display system, including direct memory 
access, will be purchased and installed in the next few months. This 
facility has been, and will be, a valuable and necessary tool for research 
into computer based control of low performance machinery at the Naval 
Postgraduate School. 

The LPC testbed is shown in Figure 2. The Controller (IBM PC- AT) is 
at the left in the figure, the Data System (IBM PC- XT) is at the right, 
and the external patch rack is in the center. The Data system can be run 
independently of the controller to generate plant modeling data, or the 
Data System and the controller can be run simultaneously through a XENIX 
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Figure 2. 



Computers For 



the Low Performance Control Lab. 
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operating system. Maximum bandwidth is about 100 Hz for single input 
single output operation. 

The Neptune manipulator is shown in Figure 3 mounted on its support 
stand. The hydraulic pump is to the right of the stand, the hydraulic 
accumulator is beneath the stand, and the bottom shelf holds the valves 
and a hydraulic manifold. The manipulator has six revolute degrees of 
freedom, operates on a water based fluid at 120 psig, and can be 
configured for solenoid ar servovalue operation. 

Besides the LPC Lab, two more or less temporary experimental 
facilities have developed cut of the NPS Robot Firefighter Project, these 
are the drilling apparatus and the flexible link facility. Key features 
of the drilling apparatus were shown in Figure 1 . The USAF drill was 
mounted on an "H" mounting and afixed to a movable milling bed platform. 
The vertical support in the "H" mounting was a steel cylinder of 1 inch 
diameter which was instrumented for strain in order to estimate the loads 
incurred by the drillnig process. The use and analysis of the drilling 
apparatus is described in detail in Section IV of this report. 

The flexible link facility is shown in Figure 4. The facility 
consists of a hydraulic pressure source (2000 psig, shown at left), a 
support platform and control box shown at the center, a one degree of 
freedom hydraulic motor actuator at right, and a flexible box beam shown 
vertically upwards. The beam was designed to flex in one plane only in 
order to study distributed dynamic structural flexibility. This facility 
was designed and implemented based upon our desire to use the firefighter 
mechanism to accomplish a range of tasks- expecially the heavy lifting 
associated with manipulating bombs, rockets, or large aircraft parts. 
Consequently, if the mechanism is not to become ridiculously bulky and 
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Figure 3. The Neptune Manipulator 
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The Flexible Link Facility* 



heavy, it can be expected to be flexible; that is, it should be designed 
to be flexible. Such a flexible mechanism poses interesting problems in 
terms of modelling and control. Much of vftat we have learned in the 
flexible link studies is being incorporated into space- based manipulator 
analysis at the time of this writing. 

B. Major Conclusions 

The conclusions which follow are referenced to the theses in which 
they were first mentioned. The theses are listed in the next portion of 
this report. 

1 . Manipulator Design for the USAF Drill 

(a) Drill walking at the inception of drilling can lead to large 
lateral forces and corresponding twisting moments being exerted on 
the manipulator support structure (Lawrence) . This can be largely 
avoided by pilot hole drilling or drill tip redesign. 

(b) Manipulator load forces vary significnatly as the drill tip 
wears (Lawrence) . Design of the manipulator must account for the 
worst case to be expected. At some point the drill becomes so 
worn that it is unable to penetrate the skin (Yobs) . A means for 
sensing the onset of this condition should, be provided in the 
system design and utilization. 

(c) The drill tip loads have been characterized for both the worn 
drill (Yobs) and the new drill (Lawrence) . Dynamic and peak loads 
have been measured. This data will allow the subsequent design of 
mechanism and acutator components for a drill manipulator. 

(d) Flexibility effects are difficult to predict and measure 

(Petroka) . This conclusion suggests, that two manipulators may be 
a good idea: one being rid id, low powered, and low performance to 
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manipulate lightweight tools such as the USAF drill; and one being 
strong, high powered, and perhaps high performance to manipulate 
heavy objects. Hie strong manipulator may be flexible and would, 
consequently, require a very sophisticated control system, perhaps 
even a vision- based controller. 

2. Rigid Manipulator Modeling 

(a) Simulation of robot motions should be done in global 
coordinates to avoid singularity problems associated with local, 
or joint coordinate formulations (McCarthy) . Such global 
simulations are easier to formulate and understand than local, 
recursive forms (Mohammed) . 

(b) For motions less than 30 degrees, the Neptune robot appears to 
be a very nearly linear machine (Harris , Lewis) . This result has 
strong implications in the design of controllers for low speed 
manipulators such as the Neptune. 

3. Rigid Manpulator Control 

(a) The USAF drill should not be used at incidence angles greater 
than about 12 degrees from the local normal, or at locations vhich 
are supported by stiffners beneath the skin (Yobs) . 

(b) Optimal Control (LQR) has been demonstrated as an appropriate 
means of achieving effective coordinated control of rigid link 
motions (Harris) . 

C. Thesis Advised 

1. McGalliard, G. R. , "A General Simulation Program for Robot 
Manipulator Arm Dynamics," NPS MSME Thesis, Sept. 1984. 

2. MaCarthy, W. F. , "Simulation of High Speed Motion of Rigid, 
Revolute Mechanisms," NPS MSME Thesis, Dec. 1985. 
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3. Lewis, D. R. , "Modelling of a Low Performance Rigid Revolute Robot 
Arm," NPS MSME Thesis, Dec. 1985. 

4. Yobs, R. L. , "Manipulator load Forces for a Robotic Firefighter," 
NPS MSME Thesis, Mar. 1986. 

5. Petroka, R. P., "Computer Simulation and Experimental Validation 
of a Dynamic Model (Equivalent Rigid Link System) on a Single-Link 
Flexible Manipulator," NPS MSME Thesis, June 1986. 

6. Burrill, L. D. , "A Feasibility Investigation for Optimal Robotic 
Control," NPS MSME Thesis, Sept. 1986. 

7. Mohammed, K. , "Non- Singular Modelling of Rigid Manipulators," NPS 
MSME Thesis, Dec. 1986. 

8. Harris, J. P., "Investigation and Development of a Micro- computer 
Based Robotic Controller," NPS MSME Thesis, June 1987. 

D. Theses in Progress 

1 . Sanders , D. , "Optimal Control as an Autonomous Vehicle Path 
Planner," NPS MSME Thesis, Dec. 1987 (expected). 

2. Fancher, C. , "Microcomputer Based Optimal Control for a Low 
Performance Manipulator," NPS MSME Thesis, Dec. 1987 (expected). 

3. Altinok, S., "Non- Singular Simulation of Rigid Manipulator. 
Motions," NPS MSME Thesis, Dec. 1987 (expected). 
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Section III 



Investigation and Development of a 
Micro- Computer Based Robotic Controller 
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situations. Deep submergence undersea exploration vehicles, unexploded 
munitions disposal units and carrier flight deck firefighting have stood out as 
examples of jobs that could be taken over by robots with the goal of reducing 
the risk to human life [Refs. 1,2,3]. In the area of shipboard fires, it was 
realized that a motorized firefighting unit cannot replace a complete 
firefighting team, but a quick reaction vehicle could provide intial fire 
suppression and quite possibly keep the fire from spreading out of control [Ref. 
1]. Additionally, the vehicle could be sent into areas where the fear of either 
ammunition or fuel explosion would prohibit a firefighting team from 
entering. 

The design and operation of a control system for such a robotic device is 
one area that has recently come under renewed scrutiny [Ref. 2], The majority 
of all robots available today use a somewhat large and immobile computer as 
the master "brain". These large computers are needed to handle voluminous 
calculations that are demanded by the complex control algorithms for 
multiple robot manipulator links. In light of their sheer size however, these 
systems are a poor choice for mobile operations. A more suitable choice is a 
micro-computer based controller. Due primarily to increases in speed and 
computational power, "personal computers” (PC’s) and the micro-processors 
that drive them, have been thrust into the limelight as the bridge to the self- 
contained robot operations. 

The first goal of controller development is the selection of the algorithm 
that will provide adequate system performance. Because cost is a major factor 
in most operations, an ideally suited controller should minimize the power 
required to conduct an operation while balancing this against the 
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I. INTRODUCTION 



Robotics, only recently considered science fiction, has now become a 
permanent part of technology. Although the term "robot” is often 
misinterpreted, robotic type devices can now be found in all aspects of society, 
from personal toys, to an integral part of many major manufacturing 
industries. For the purpose of this paper, a robot is defined as a self contained 
electro-mechanical system consisting of: a digital controller, a power source, 
and a manipulator arm. 

I 

The advantages of using machines to replace man /'have long been 
appreciated. In industry, economics has driven the incorporation of this new 
technology to accomodate and sometimes replace older methods. The 
resulting automation has markedly boosted output while decreasing the per 
unit manpower costs. This has in turn caused an overall increase in 
production efficiency making many items more affordable to the consumer. 
Another advantage of automation has been the attainment of more consistent 
results. Painted items thus become more uniform, machined items become 
more exact, and production quality becomes more standard from item to item. 
All of these facts not withstanding, however, the most important reason for 
the use of mechanical equipment is to lessen the endangerment to man while 
operating in hazardous environments. The more autonomous the robotic 
devices, the more removed and hence safer the person controlling the 
mechanism. 

Recently, the U.S. Navy has developed an interest in the use of robots 
with an aim toward replacing humans in potentially life threatening 
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positional error generated by the control system. To accomplish this, a 
modern controller should be investigated for implementation. Such modern 
controllers are capable of coordinating multiple joint inputs to control the 
complicated manipulator output configuration. The proof of these controllers 
lies in their ability to control a real piece of hardware. 

This thesis describes the implementation of a modern controller on the 
Naval Postgraduate School manipulator control test bed. The test bed consists 
of: a NEPTUNE multiple rigid link robot system (Fig. 1), an IBM PC-AT 
micro-computer, and two Atchley jet pipe electro-hydraulic servovalves. The 
NEPTUNE system was modified so that it could be run as originally 
purchased (solenoid operation), or through a PC based servo controller in 
order to demonstrate the comparative performance associated with each type 
of controller. 
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Figure 1 

NEPTUNE Manipulator 
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II. BACKGRQ UND/PROBLEM DESCRIPTION 



The primary goal of this work was to initiate an NPS test bed system for 
Navy development of robotic controllers. The fact that operational robots will 
be needed to operate in a hostile environment suggests that water based 
actuators be used as the prime motive force. The NEPTUNE robot arm and 
control system fit this requirement and were previously selected [Ref. 4]as the 
system nucleus. Lewis [Ref.4] described the operation of the basic 
computer/solenoid controlled assembly which was modified in the present 
investigation. 

In the original NEPTUNE actuation scheme, a pseudo-proportional 
control was obtained through the use of several fluid restrictors (Fig. 2). By 
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Figure 2 

Solenoid hydraulic system 
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ganging more than one of these restrictors together on a joint, the motion and 
flow rate history as seen in Figure 3 was acheived. However, in order to 
investigate advanced servo controls, an infinitely variable flow rate is 
required as in Figure 4. In order to modify the NEPTUNE to accept a servo 




Figure 3 

Solenoid operation, (a) motion, (b) flowrate history. 




Figure 4 

Servo valve operation, (a) motion, (b) flowrate history 
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control signal, Burrill [Ref. 5] conducted a feasibility study to determine the 
plant requirements and select the best suited servo valve. "Off the shelf' 
valves were desired to both minimize the initial procurement costs, and also to 
standardize the system components with those generally available to 
industry. The Atchley model 218 jet pipe servo valve was selected, and two 
valves were placed in service. Because only two dimensional control was to be 
demonstrated, axes 1 and 2 were chosen (Fig. 5). 



An initial hardware modification was required to improve hydraulic 
system operation. To provide a reservoir of fluid energy, the NEPTUNE 
hydraulic power pack comes with a 42 cubic inch air bladder accumulator. In 
order to reduce entrained air in the water-based operating fluid, an air 
bladder in the accumulator is filled to approximately 55 psi, or one half of 
system operating pressure. The fluid pump then fills the accumulator and 
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Figure 5 

NEPTUNE axes of rotation 
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shuts off at 110 psi. As fluid is drawn down in an arm move, an unloading 
valve monitors the pressure and automatically turns the pump on at 
approximately 100 psi. With two joints being operated simultaneously, pump 
pressure waves can be observed in the actuator cylinder data. This effect was 
difficult to predict and evaluate. 

To alleviate this problem, joints one and two cylinder dimensional 
measurements were taken, and the maximum amount of fluid needed to move 
both actuators was calculated. It was determined that approximately 36 cubic 
inches of fluid was needed to move both actuators the length of their travel. A 
580 cubic inch accumulator was fitted such that less than 10 percent of the 
total volume would be used on any arbitrary move. The resulting cylinder 
pressure data was much more well behaved. 

At this point some clarifications should be made as to the distinction 
between the plant and controller as described by this research. Figure 6 
shows the general multi-port layout with the term system defined as the 
combination of the computer, both servovalves, and their associated links. 
The term plant refers to the arm and servovalve combination. Although not a 
radical departure from conventional thinking, the inclusion of the servovalve 
in the plant block instead of the controller block could be considered unusual. 
There were two reasons for taking this approach. First, the plant is normally 
considered as that part of the system that is, for one reason or another, already 
selected, unchangeable and must be controlled. Secondly, as will be discussed 
later, the dynamic characteristics of the links when combined with the 
hydraulic characteristics of the valves allows a simple yet effective modeling 
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Ps = constant — ► SUPPLY 
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Figure 6 

Multi- Port Block Diagram 




of the plant. Because the speed and computational efficiency of the IBM PC 
AT are significantly less than that of a larger mainframe computer, any 
simplification of the model may allow for a shorter controller sample rate.The 
controller is defined as the combination the sensor acquisition/conversion 
elements, the IBM PC AT computer, and the control program. With the 
exception of Digital-to-Analog Converters (DAC), Analog-to-Digital 
Converters (ADC), and the computer itself, all control and regulation is 
software based. The DAC and ADC are contained on two expansion port 
boards internal to the AT. Computer interfacing was done with the system 
known as LABPAC, which has 6 12 bit DAC's and 16 12 bit ADC's. With 12 
bits of data availabe for conversion, the maximum available resolution is .024 
%. In addition to the converters, LABPAC has an onboard clock timer for 
sample rate regulation. Once setup, LABPAC automatically does both data 
acquisition and control at set intervals. 

The system model is based on the Multiple Input Multiple Output 
(MIMO) elements of: the computer controller, servovalves, hydraulic 
actuators, and the arm itself. To create motion in the system, a voltage 
potential (v) is applied to the servovalve electrical input, which in turn sets up 
a current flow through the valve actuating coil. This causes an opening of the 
valve orifice. The actuating voltage in combination with the actuator 
feedback of fluid flowrate (q) produces the valve pressure outputs Pi and P 2 . 
In turn Pi and P 2 are used as inputs to the actuator along with the arm 
velocity (0) in order to generate the torque (T) which creates motion in the 
arm. 
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In order to validate the robot test bed system for future operation, an 
optimal state space controller (LQR) was developed and implemented. As 
discussed by Ogata [Ref. 6] and Owens [Ref. 7], optimal control theory 
requires that a merit function J be selected. J can be considered the 
performance index that is to be minimized. To meet with established 
conventions, the merit function usually uses the two quadratic forms that 
follow. 



</j= (x T Qx) dt ( Eqn . 1) 

J 2 = (u 7 Ru) dt (Eqn. 2) 

The matrices x and u are the state space matrices that correspond to state 
error and system-inputs respectively. Q is a weighting matrix that determines 
the relative importance of each term in the state matrix. R is determined by 
the system performance factor to be minimized, and cost weighting given each 
input. J can take the form of the arithmetic combination of Equations 1 and 2 
in the form : 

•/ = •/. (Eqn. 3) 

J = Jq {A (x r Qx) + (u r Ru)}cft (Eqn. 4) 

This is helpful because for most real life systems the minimization of the 
state error alone must be tempered by at least two considerations. First, all 
operations "cost" something whether it is energy, time, or some other 
quantity. Secondly, since this form of optimal control produces the feedback 
law: 
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u = -Kx 



(Eqn. 5) 



a large state error will consequently produce a correspondingly large 
correcting input. Too large a feedback, however, can saturate most real 
systems thereby giving them a high degree of non-linearly. 

To narrow the scope of the problem and provide a realistic goal, the 
performance indices chosen were related to state error versus the square of the 
hydraulic power used by the system. For the purpose of this discussion, the 
definition of hydraulic power will be that used by Merritt [Ref. 8] where: 

Power = q (Eqn. 6) 

P = P , — constant 

s supply 

Initial dynamic modeling of the two link arm mechanism was done in 
global vice local coordinate systems. This decision was made in anticipation of 
future problems with solving for joint accelerations using the Newton-Euler 
forward dynamic equations. If local coordinates are used in conjunction with 
the Denavit-Hartenburg [Ref. 9] convention, the potential exists for a 
singular solution when two adjoining links are colinear with the relative 
angle being either 0 or 180. degrees It was realized that the joint solutions in 
local coordinate systems are more computationally efficient; However, the 
singularity condition described above precludes general implementation of a 
model reference position control algorithm that does not have error handling 
procedures. 



27 



Subsequent modeling of the NEPTUNE system revealed that a much less 
complex plant model was needed to satisfactorily describe the motion. This 
realization, when coupled with the fact that all sensors and actuators 
functioned in local coordinates, allowed final implementation of the controller 
in joint relative coordinates. 

There have been numerous research efforts aimed at various types of 
micro-computer robotic control [Refs. 10-14]. Optimal path control is dealt 
with by Snyder [Ref. 15]. The closest research done in the area of optimizing 
power appears in Reference 16, where performance indices similar to the one 
mentioned above are used, with the focus on weighted near minimum time- 
fuel control. 
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m. ANALYTICAL SYSTEM CHARACTERIZATION 



A. RIGID TWO LINK SYSTEM MODEL 

The Newton-Euler formulation for the robot arm was used in developing 

the classical dynamic model. Previous work by Lewis [Ref. 4] on the 
NEPTUNE arm showed that with the original solenoid valves installed, 
reduced dynamic Equations 7 and 8 could be used (all angles are relative). 

T l = © 2 + Dg sin(& 2 — ©^ + G g siniOJ — M ^ ( Eqn . 7) 

T 2 = P 2 © 2 + Dgsin(& 2 — © t ) (Eqn. 8) 

D = m 2 d 2 + m 3 (l 2 + d 3 ) 

G = (m. } m.j) + m, d^ 

\f = return spring moment 

. S 

Because the modified arm system contained two servo valves of a higher 
performance and flow rate, the full Newton-Euler two degree of freedom 
equations of motion were used [Ref. 17]. The global formulation is shown in 
Figure 7. 

The torque required for each joint was solved for as follows : 

7\ = (J { + cf 2 + m 2 1* )Q l + m 2 l { d 2 cos <© t - ©^C^ 

+ m 2 l l d 2 sin(Q l -Q^& 2 +(m. 2 g + m y g d ] )cosQ l + T 2 (Eqn. 13) 
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r 2 = (J 2 + m 2 dl)Q 2 +m 2 l l d 2 cos(e 2 -Q l )Q l 

+ m 2 l l d 2 sin(Q 2 -Q l )Q* +m 2 gd 2 cos& 2 (Eqn. 14) 




B. SERVOVALVE MODEL 

The servo valve was modeled using the classical equation for flowrate: 

q = k iV P * (Eqn. 15) 

v 

1 ~ Q 

Since all inputs to the actuating coil were DC, the valve impedance was 
modeled as a simple coil resistance (Q). - This allowed a direct conversion from 
the model, in current (i) to actuating voltage (v) (Eqn. 16). With the servo 
valve installed in the system, the flowrate relationships in Figure 8 could be 
used determine either the flow or pressure drops at any point in the system, if 
other values are known. 
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Figure 8 

Pressure/Flow Schematic 

If the valve is assumed to be symmetrical (both ports and spool areas are 
identical), then the valve constant k will be the same for both supply and 
return. Once k is determined experimentally , Equation 16 can be rearranged 
as follows : 




{Eqn. 17 ) 



Refering to Figure 8 and Figure 2, it can be shown that 0 and the flowrate q 
are linearly related by: 



q. = C0 


(Eqn. 18) 


in 3 

c = .1 131 — 

deg 


(Eqn. 19) 




because of the combination piston, and rack and pinion gear assembly. With 
P 3 = constant, Equation 17 can be further reduced to: 



C. VALVE TESTING 

Pressure/flowrate tests were conducted on both servovalves in order to 
determine the valve constant k (Table I). As can clearly be seen, although k 
did follow a trend for each run, k was generally not a constant-either between 
the ports of the same valve, or with different input voltages. One explanation 
for this is that the normal operating pressure range for this valve is between 
200-3000 psi. The pump and accumulator currently operate at between 100- 
110 psi. It is possible that the present configuration of the system forces the 
valve to operate in a non-linear region below the minimum design point. This 
notion is bolstered by the fact that although the rated maximum input voltage 
is 8.0 VDC, the observed saturation point is approximately 1.0 VDC for joint 
1,- and 0.6 VDC for joint 2. At inputs greater than these no corresponding 
increase in flowrate occurs. 

Pressure plots of the driving cylinder (P t ), and the opposing cylinder (P 2 ) 
showed expected trends for a 0 to 90 degree move of link 2 (Fig. 9). For link 1, 
Pi and P 2 remain at equal pressures until ©i became greater than 35 degrees 
(Fig. 10). This is believed to be caused by the driving force on link 1 of the 
return springs. Once ©i exceeds 35 degrees, the force added by the springs 




( Eqn . 20 ) 



or 



P = fen (v , 0 ) 



(Eqn. 21 ) 
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TABLE I 



VALVE CONSTANT CALCULATION RESULTS 





Input 


Flow 


delta 


delta 


Ks 




current 


rate 


Ps 


Pr 


- 




(mA) 


(in/s) 


(psi) 


(psi) 




Joint 1 


2.5 


1.43 


65.9 


64.7 


70.4 


Run 1 




1.15 


54.0 


74.1 


62.4 






.855 


42.0 


87.6 


52.8 


Joint 1 


5.0 


1.86 


67.4 


65.0 


45.3 


Run 2 




1.47 


56.7 


75.4 


39.0 






1.07 


43.6 


90.3 


32.4 


Joint 2 


1.25 


0.410 


40.7 


101.3 


51.4 


Run 1 




0.446 


47.4 


98.2 


51.8 






0.541 


61.0 


92.7 


55.4 


Joint 2 


3.75 


0.95 


36.0 


96.3 


42.2 


Run 2 




1.12 


40.8 


90.1 


46.8 






1.36 


49.5 


82.5 


51.6 



Kr 



70.4 
53.3 

36.5 

46.1 

33.9 

22.5 

32.6 

36.0 

45.0 

25.8 
31.5 

39.9 
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Figure 9 

Joint 1 Pressure Plot 
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THETA I (DEG) 




Figure 10 

Joint 2 PressurePlot 
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THETA 2 (DEG) 




decreased significantly, and the pressure plots follow the anticipated pattern. 
The maximum change in either or P 2 , for either link was determined to be 
25 psi. When incorporated into the valve model (Eqn. 15), the change in the 
flowrate or joint velocity is relatively small, even with a large angular, 
change. These observations tend to confirm Lewis’ hypothesis that the 
NEPTUNE is a low performance robot arm, and as such, most of the system 
dynamics may be discounted. 

Because classical valve operation and arm dynamics were unable to be 
validated by actual component performance, it was decided that an analytical 
model of the system would provide poor results. As an alternative, an 
equivalent linearization technique utilizing frequency response was selected. 
The system was modeled empirically using magnitude and phase plot. 
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IV. EMPIRICAL SYSTEM CHARACTERIZATION 



A. EXPERIMENTAL SETUP 

A Hewlett-Packard 3562A Dynamic Signal Analyzer was used to produce 
the Bode plots for plant modeling and analysis. This signal analyzer is 
capable of scanning a desired frequency range with a sinusoidal output, while 
simultaneously measuring and plotting transfer function magnitude and 
phase data. The signal analyzer output was hooked directly to the servovalve 
electrical inputs. The input to the signal analyzer was taken from the existing 
position feedback potentiometers on the manipulator arm. A small amplitude 
voltage was used to excite the servovalves. In this way an equivalent 
linearization of the arm was measured. 

Based on the most likely operating regions of the arm, three linearization 
configurations were selected. For link 1, a mid-range point of 45 degrees was 
chosen. It was believed this position best accounted for actuator dynamics, 
link motion, and the nonlinearity created by the two return springs on link 1. 
Three linearization points were chosen for link 2. They were -45, 0, and 4-45 
degrees (Fig. 11). These positions provided the best coverage for link 2 
motion. Table II lists all data acquisition runs and initial conditions. Figure 
12 displays the anticipated MIMO plant model to include cross coupled 
dynamics between the joints. It was anticipated that each linearization point 
would require a different set of transfer functions. 

Although the effects of servovalve fluid leakage between the spool and 
housing were initially deemed negligible compared to the flow through the 
valve, the experiments showed that drift in the position of the link became 
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TABLE II 

EQUIVALENT LINEARIZATION DATA RUNS 



Run No. 


Ui 


u 2 


01 


©2 


Function 




(VDC/DC 


(VDC/DC 


(deg) 


(deg) 


Sampled 




offset) 


offset) 


1 


.3/.1 


fixed 


+ 45 


-45 


© i/Ui 


2 


.3/.1 


fixed 


+ 45 


0 


©i/Ui 


3 


.“H.l 


fixed 


+ 45 


+ 45 


©i/Ui 


4 


fixed 


.15/. 1 


+ 45 


-45 


©2/U 2 


5 


fixed 


.4/. 15 


+ 45 


0 


0 2/U‘2 


6 


fixed 


.2/.1 


+ 45 


+ 45 


© 2 /U 2 


7 


fixed 


.15/. 1 


+ 45 


-45 


©l/u 2 


8 


fixed 


.4/.15 


+ 45 


0 


©i/u 2 


9 


fixed 


.2/.1 


+ 45 


+ 45 


©i/u 2 


10 


.3/.1 


fixed 


+ 45 


-45 


0 2/Ui 


11 


.3/.1 


fixed 


+ 45 


0 


02/Ui 


12 


.21.1 


fixed 


+ 45 


+ 45 


02/Ui 



Ol= 45’ 
©2 = -45’ 



©1= 45’ 
© 2 = 0 



©1=45 
©2 = 45 




Figure 11 

Linearization configurations 
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Figure 12 

Plant Block Diagram 

noticeable when the valve was hydraulically charged, but electrically 
inactive. To negate this effect, the unactivated link fluid supply hoses were 
disconnected and plugged at the valve end. The result of this was that the 
hydraulic characteristics (bulk modulus of elasticity, fluid properties, 
pressure wave reverberations, etc.) of the actuating cylinders and supply 
hoses on the unactivated link was included in the equivalent linearization, 
with the servovalve assumed to have no leakage. 
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B. RESULTS 



1. Link 1 

Figure 13 shows the results of the frequency response plots at 0i = 
+ 45 degrees, 02= -45 degrees. In this run, 02 was held constant and the 
sinusoid signal was input to servovalve 1. The output from the joint 1 feedback 
potentiometer was compared to the valve input signal in order to obtain the 
Bode plots. The frequency range of measured output was from .1-10 Hz. Curve 
approximations to data (dashed line) showed the best fit in both magnitude 
and phase plots was a two pole system. One pole was located at 0 Hz and the 
second pole was rooted at 3.6 Hz. The gain K p was 1.1. The system displayed 
an apparent transportation lag of 0.04 seconds. This feature was anticipated 
since all real systems, particularly mechanical ones, have some finite time 
delay between signal input and noticeable output. The transportation lag 
was constant for all data runs of both valves, mainly because both valve and 
cylinder assemblies are essentially identical with the supply hoses of about 
the same length. A trial first order curve fit of the same data with a single 
pole located at 0 Hz corresponding to a transfer function of (K/s) was also 
tried. Although the approximation was relatively good at lower frequencies, 
higher frequencies showed poor correlation in both magnitude and phase. 
Adjustments for transportation lag provided no significant improvement and 
the overall result was considered unacceptable. After careful analysis, the 
following transfer function was determined to most closely model link 1: 

1 . 1 e ~ 004 4 



( s + 3.6) 



(Eqn. 22) 



0i/Ui 




lOO. Q2m Hz ID 



< 




Figure 13 

Configuration 1 Linearization - 0i/Ui data 
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Linearization of link 1 at 02 = 0 degrees (Fig. 14), and + 45 degrees (Fig. 
15) produced similar results as described above with the transfer functions 
being : 



0 i l.le -0 04a 
C/j ~ s( s + 3.6) 



(Eqn. 23) 



0 ! _ 1.0 e -004s 
U x ~ s(s + 3.6) 



(Eqn. 24) 



respectively. Since all three transfer functions were within an experimental 
errror factor of 10%, it was decided that one common linearization point would 
suffice (Eqn. 22). Converting 0i from volts/sec to degrees/sec gave the final 
form that would be implemented for the link 1 transfer function. 



0 ! 100.16 e ~ 0 - 04 * 

U ^ s ( $ 4 - 327 . 7) 



{Eqn. 25) 



2. Link 2 

Figure 16 shows data acquired from link 2 with ©i fixed at +45 
degrees, and ©2 linearized around -45 degrees. Analysis of the plots revealed a 
pole at 0 Hz, similar to link 1 results. However, no signs of any higher 
frequency poles were evident. A good curve fit to data was accomplished with 
a transfer function of the form (K/s). Once again a transportation lag of 0.04 
seconds was evident in the phase plot, and the gain K p was determined to be 
0 . 21 : 
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Figure 14 

Configuration 2 Linearization - 0i/U i data 
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0 1/Ui 





lOO. Q2m Hz lO 



Figure 15 

Configuration 3 Linearization - 0i/U i data 
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Figure 16 

Configuration 1 Linearization - 02 /U 2 data 
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e 2 .21 e -0 04 * 

— = (Eqn. 26) 

U 2 S 

Figure 17 shows the results of linearization around 02 = 0 degrees. 
Comparable results were obtained from this data run, with the transfer 
function determined to be : 

e 2 .26 e _004s 

— = (Eqn. 27) 

U £ s 

Linearization at 02 = +45 degrees could not be accomplished using the 
signal analyzer because of the low DC signal offset required to keep the link 
in gravitational equilibrium. Limited data points were taken manually with 
Fig. 18 showing an approximate curve fit. Although this data was unusable in 
providing an accurate model, it does show correlation to the data of the 
previous two linearization points. The approximate transfer function for link 
2 was: 



0 

U. 



.23 e 



-0.04 s 



S 



{Eqn. 28) 



As in the case of link 1, since the results were so closely grouped, it was 
felt that only one linearization equation would be needed to adequately model 
the range of the three original set points. Conversion from volts/sec to degree/ 
sec produced the final form for the link 2 transfer function (Eqn. 29). 
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Figure 17 

Configuration 2 Linearization - 02 /U 2 data 
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(Eqn. 29) 



0 2 {22.5 e~ 0 04 5 } 




3. Cross Joint Coupling 

In most complex mechanisms, dynamic cross coupling between links 
is a very real factor that will inevitably increase the system complexity. In an 
effort to account for this the plant model presented in Figure 12 contains two 
additional transfer functions, 0i/U2 and © 2 /U 1 . These transfer functions 
contain the cross-coupling terms that generate an output motion in link 2 
when joint 1 is actuated, and vice versa. Data taken in runs 7 through 12 
were used to determine the cross-coupled effects. In all runs, random data 
scatter occured between -35 and -70 dB with no correlation or consistency. The 
conclusion is that the low velocities and accelerations generated by this arm 
are not enough to produce any significant cross-coupling motion. In effect each 
individual link can be considered independent. This result is specific to the 
NEPTUNE arm performance and may not be the case for higher performance 
systems. Therefore, care must be taken to test for the existence of dynamic 
cross-coupling before attempting similar simplification. 

4. Revised Plant Model 

With the elimination of cross-coupling and the decrease in the 
number of linearization points, a revised plant model is shown in Figure 19. 
The conversion of these transfer functions from the Laplace domain to the 
time domain yielded : 



= —327.7 + 100.2 U^it— 0.04) (Eqn. 30) 

0 2 = 22.5 U 2 (t- 0.04) (Eqrc.31) 
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The acceleration term in Equation 30 strongly suggests that although the U 2 
input has no cross-coupling effect to © 1 , the inertia generated by the mass 
translation of link 2 is still important to the motion of link 1. The impact of 
Equation 31 is that the velocity of link 2 is directly controllable by the input 
voltage U 2 , to the servovalve. 




Figure 19 

Revised Plant Model 
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V. OPTIMAL CONTROL THEORY 



A. CONTINUOUS TIME LINEAR QUADRATIC REGULATOR 

The Linear Quadratic Regulator (LQR) is considered one of the most 
important and most used optimal controllers. The reason for this is based on 
the fact that although the derivation of the control law is a laborious and 
painstaking procedure, the implementation is an extremely simple and 
straightforward method to achieve multiple coordinated inputs. LQR control 
assumes that the plant can be modeled linearly in state space with: 

x = Ax + Bu (Eqn. 32) 

Additionally, the performance index J must be of the quadratic form: 

J = /*{.\ (x T Qx) (\x T Ru)} dt iEqn 4) 

As discussed earlier, the weighting matrices Q and R have special 
definitions with respect to state error and control cost. Further restrictions 
apply to each weighting matrix : Q must be either a positive semi-definite or 
positive definite, real symmetric matrix; R on the other hand is required to be 
a positive definite, real symmetric matrix. Once Q and R are reduced to the 
proper form, a gain matrix K can be found in a two step procedure. The first 
step necessiates the solving of the reduced matrix Riccati equation for the 
unique positive definite, real symmetric matrix P. 

A T P + PA-PBR' 1 B T I , + Q = n iEqn. 3 3) 

The matrix K will then follow directly as: 
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K = R~ l B r P 



( Eqn . 34) 



With K determined, the governing control law is given as 



u = — Kx iEqn. 5) 

Figure 20 shows the corresponding system block diagram. 




Figure 20 

Optimal Feedback Control Block Diagram 



In this work the state selections were: xi = Oi X 2 = ©i t X 3 = 02. Since 

# 

NEPTUNE has only position feedback potentiometers, joint 1 velocities (0i) 
were calculated using a backward difference method. 

B. DISCRETE TIME LINEAR QUADRATIC REGULATOR 

Because the system controller was to be implemented in software using a 
computer, the continuous time LQR had to be converted to a discrete time 
format. This required discretizing the plant model, the performance index, 
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and the matrix Riccati equation. This work has been done previously [Refs. 
4,15,16] and the results are presented in Equations 35 (plant), 38 
(performance index), and 39 (gain matrix) for comparison to the continuous 
time LQR controller equations. 



x(n +• l ) = <t> x <n) + I' u (m) (Eqn. 35) 

4> = e A/l (Eqn. 36) 

r = {/q e A * dt } B (Eqn. 37) 

•J =^> {.\ (x T (n)Qx(n ) )+ (u^Ca) Ru)(n)}d< (Eqn. 38) 

K = (R+r T Prr‘r r PD (Eqn. 39) 



C. EFFECTS OF TRANSPORTATION LAG 

Transportation lag in the discrete time system is an important factor from 

the viewpoint of sample rate and system stability. If the sample period Ih) is 
significantly shorter than the time delay (x), several control periods may go by 
before there is a corresponding change in the output. One method to rectify 
this problem is provided for in Reference 18; It requires the addition of an 
extra state for each overlapped sample period. With t less than the sample 
period, Equations 40,41,42,43 detail the procedure and contents of the 
modified $ and r matricies. Actual determination of the sample period that 
was used will be discussed in a later section. 
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0 = e 



(Eqn. 40-) 



\h 



r x = e Alh -" {/Jc At dt } B {Eqn. 41) 

r 0 ={/' 0 A - x ’ e A ‘ dt) B {Eqn. 42) 



x (n+ 1 ) 


,, 


<p 


x(n) 


4- 


r o 


u < n) 




0 0 


. U (ft— 1) 




i 



{Eqn. 43) 



D. DISCRETE PLANT DETERMINATION 

From Equations 30,31 the continuous time linear plant model was 

determined to be: 



— 












i 


0 


0 




0 


0 


0 


-327.7. 


0 


X + 


100.2 


0 


0 


0 


0 




0 • 


22.5 














_i 



(Eqn. 44) 



Discretization of this for a sample period of 0.05 seconds leads to: 



x (n-f- 1) 



1 


3.05(10-3) 


0 




1.44(10-2) 0 


0 


7.66(10-8) 


0 


X(n) + 


3.03(10-1) o 


0 


0 


1 




0 1.125 



(Eqn. 45) 



Inclusion of the transportation delay produced the following $ and V matrices: 



* = 



1 

0 

0 

0 

0 



3.05(10-3) 

7.66(10-8) 

0 

0 

0 



2.03(10-2) 

1.15(10-2) 

0 

0 

0 



0 

0 

0.9 

0 

0 



( Eqn. 46) 
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r 



( Eqn 47) 



2.16(10-3) 

2.94(10*‘> 

0 

1 

0 



0 

0 

2.25(10-1) 

0 

l 



E. DETERMINATION OF THE WEIGHTING MATRICES 

The selection of the R matrix was driven by the goal to minimize 

hydraulic power. To fit the form of the quadratic performance index and also 
to force R to be positive definite it was decided to minimize Power 2 . 



Power ~ = ip q.) 

■i *> 



(Eqn. 48) 



. = P* c 1 ©? 

J # 



= Ru? 



Joint 2 fit the form exactly with 02= 22.5 U-2- Because of the second order 
term in the linearized equation for joint 1, however, a similar form did not 
exist for 0i andUi. To alleviate this, a reversion back to the first order 
approximation ofjoint 1 motion was required for the purpose of merit function 
evaluation only . As discussed earlier, a satisfactory fit at lower frequencies 
was found for the (K,s) joint I model, but increasing deviation at higher 
frequencies was observed. Actual observation of the link in operation 
indicated that the higher frequency dynamics occured mainly in the 
transition from static equilibrium to sustained motion, and vice versa. For the 
majority of the operating region much lower velocities were encountered. 
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Based on this it was felt that, although the approximation was not accurate 
over the entire operating spectrum, reasonable results could be obtained for 
power comparisons with the (K/s) model for link 1. This led to the following R 
matrix based on first order models for both joints: 



80473 0 

0 78400 



(Eqn. 49) 



The Q state weighting matrix took the form of a reduced identity matrix 
I. The states that were considered important were Oi, Oi, and 02- Minimizing 
the feedback states of past inputs, or Ui(n-l), was felt to have no real utility 



Q = 



1 0 0 0 0 
0 10 0 0 
0 0 10 0 
0 0 0 0 0 
0 0 0 0 0 



(Eqn. 50) 



and hence was ignored. Although discussed only briefly in most control texts, 
the proper selection of Q is the major factor in determining which states the 
controller will emphasize in optimizing plant performance. 

In order to ascertain the affect that various magnitudes of X (Eqn. 38) 
have on power minimization, three trial values were chosen. The first value of 
30, was picked specifically to show a situation where both valve 1 and 2 intial 
voltages were in the saturated region. The second value of 15 placed the intial 
voltages at the saturation boundary. The final value of 1 produced valve 
inputs that remained in the acceptable performance range throughout the 
move. 
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F. CALCULATION OF THE GAIN MATRIX K 

With <&, T, R, and Q determined, the Riccati equation was solved 
The following values were obtained. 





3.52(1 0 _3 ) 


1 08(10“ 5 ) 


0 


7.17(10"°) 


0 


K .v=l = 


0 


0 


3.56(1 0 _ 3 ) 


0 


3.2(10 ~ 3 ) 




1.37(10 ~ 2 ) 


4.17(10 ~ 5 ) 


0 


2.78(10 ~ 4 ) 


0 


K ^15 = 


0 


0 


1 .37(10 _2 ) 


0 


1 .24(1 0 ~ 2 ) 




1.93(10 ~ 2 ) 


5.89(10 — 5 ) 


0 


3.94(10 " 4 ) 


0 


K V--30 = 


0 


0 


1.93(10 ~ 2 ) 


0 


1.74(10“ 2 ) 



for K. 



( Eqn . 5 1 ) 



(Eqn. 52) 



(Eqn. 53) 
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VI. MICROCOMPUTER IMPLEMENTATION 



A. CONTROLLER SAMPLE RATE 

Determination of the correct controller sample rate for the IBM PC AT in 
conjunction with the NEPTUNE arm was a two part process. First, an offline 
software test was conducted to determine the time required for the execution 
of the complete control algorithm. After this was accomplished, an 
examination of the arm model was performed to predict the minimum sample 
rate that was needed perform the ADC, control, and DAC while still 
maintaining system stability. 

The offline software test consisted of coding the entire control loop in 
Fortran, and compiling. This test code was then iterated 5000 times, and the 
elapsed time was recorded. The results of this test produced an average 
computer control cycle period of 2.57 milli-seconds. As a general rule of 
thumb, it was decided that the control computation time should not exceed 
60% of the computer control sample period. The reasoning behind this was 
that system overhead time would be required for sensor acquisition, 
conversion, and computer "housekeeping" chores. For LABPAC this overhead 
time is approximately 0.8 milli-seconds, making the minimum computer 
sample period 3.37 milli-seconds. 

Shannon sampling criteria was used to determine the minimum sample 
frequency required to enable reconstruction of the sampled motion [Ref. 11. 
Since the highest frequency pole was associated with joint 1 at 3.6 Hz, the 
lowest acceptable sample frequency for the system was set at 7.2 Hz, with a 
recommended sample frequency of 36 Hz. This equated to a sample period for 
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system stability of between .028 (recommended) and .139 (required) seconds. 
Use of the recommended sample rate, however, meant that the transportation 
lag of 0.04 seconds would extend over two sample periods. A compromise .05 
second (50 msec) sample period (20 Hz) was implemented, thus making the 
duty cycle (3.37 msec) less than 10% of the sample rate. 

B. POWER DATA COMPARISONS 

A total of four comparison runs were performed with the arm. The first 
three runs used the gain matricies that corresponded to the three previously 
discussed values of X. The last run utilized the solenoid control as originally 
configured on the NEPTUNE to provide a baseline for comparison. The start 
point for each run was 0i = 2O degrees, ©2 = -60 degrees; with the end point 
being 0i =70 degrees, ©2 = 30 degrees. Figures 21,22 show time-motion plots 
of joints 1 and 2. Each link shows expected results, with the optimal control 
runs approaching their desired set points assymtotically. It should be noted 
here that in the case of X =1, the target position was never actually reached. 
©2 finally ended 3 degrees short of the desired position, because valve input - 
Kx was in the deadband voltage for that small error condition . In 
comparison, the solenoid runs approached the target positions with constant 
velocity, and motion was slowed by use of the fluid restrictors. When the link 
was within the error deadband, motion was stopped by shutting the solenoid 
valve. 

Figures 23,24,25 compare total hydraulic power consumption for the arm 
( S(Pi l -P2 i )q i ) with solenoid power consuption, for the specified A's. The 
negative power values that appear for the solenoid and \=1 move, are 
attributable to the fact that the joint 1 return springs were providing all forces 
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Figure 21 

Joint 1 position plot 
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Figure 22 

Joint 2 position plot 
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Figure 23 

Total arm power - X = 1 
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Figure 24 

Total arm power - A= 15 



62 







CJ 

— 

cn 



£2 



0*001 0*09 0*09 0*0* 0*02 0*0 0*02- o*o*- 

(oas/Ni-ai) aajuod 



Figure 25 

Total arm power - \ = 30 
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necessary for motion of link 1 until 0i exceeded 35 degrees. This phenomenon 
was not evident for A = 15, or A = 30. As A increased from 1 to 30, the 
maximum instantaneous power used by the arm was correspondingly higher 
reflecting increased weight on error. Accurate values for the actual power 
used by the arm were difficult to determine, largely because of the noise from 
the position feedback potentiometers, and random errors in ADC. However, a 
general trend can be discerned that shows A = 1 had a lower maximum power 
requirement than did the solenoid. The sawtooth like curve of the solenoid 
operation correlates with fluid restrictor closing. The power requirements for 
all optimal control moves appear to transition smoother, from high to low 
values. This can be advantageous, especially in hydraulic systems, where 
sudden shocks are not desired. 

Total system power plots (power taken from accumulator with Ps 
assumed constant) show results consistent with arm power results . A = 30 and 
15 plots (Fig. 26,27) show that the power used by the system always exceeds 
solenoid power requirements. In contrast, the power requirements for A=1 
somewhat closely matches that of the solenoid (Fig. 28). The major deviation 
is the initial power surge for the servovalve operation. Once again, 
transitions in power levels were smoother in optimal control than solenoid 
control. 
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Figure 26 

Total system power - \ = 30 
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Figure 27 

Total system power - A. = 15 
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Figure 28 

Total system power - \ = 1 
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vn. CONCLUSIONS 



A. SYSTEM LINEARITY 

The intial evaluation of the equations of motion for the NEPTUNE 
system were based on the assumption that the rotational inertia of the arm 
structure was dominant. This in fact was not the case with the NEPTUNE. 
With the servovalves installed, the dominant mechanism for determining 
motion, was the combination servovalve and actuator assembly. Even at 
relatively low pressures, the NEPTUNE behaves linearly. In essence, the 
hydraulic system has enough power to override much of the low performance 
dynamics that might have been expected. 

B. SERVOVALVE PERFORMANCE 

The data that was obtained for the servovalve did not reflect the 
anticipated classical valve model. The reason for this is believed to be the 
relatively low system pressure and the distance between the servovalves and 
pressure sensors. It was felt that by correcting these factors, better 
correlation may be obtained. Nevertheless, equivalent linearization 
techniques were used to adequately model both the valve and arm 
performance as configured in the present study. 

C. MICRO-COMPUTER OPTIMAL CONTROL 

Micro-computer optimal control is a simple to implement, feasible control 
scheme. The control system can be "tuned" to the desired performance level 
by the proper choice of weighting matrices Q and R. Additionally, a smooth 
transition in motion and hydraulic power requirements can be expected. 
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D. SAMPLE RATE 



Since the NEPTUNE was of such low performance, the maximum sample 
rate available with the IBM PC-AT was adequate. The driving concern in 
sample rate determination was the proper handling of the transportation 
delay. Finally, although the selected sample rate was less than a 
recommended five times the Shannon sampling criteria, system stability was 
not affected. 
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VIIL RECOMMENDATIONS 



Continue the use of NEPTUNE plus IBM PC-AT as a test bed for 
controller development. Three future projects are recommended as follow on 
work. First, a three dimensional optimal controller should be implemented in 
order to utilize the full range of motion on the NEPTUNE arm. Second, 
comparison studies of optimal versus Proportional-Integral-Derivative (PID) 
control should be conducted, to provide a better baseline evaluation against a 
common type of controller found in industry. Finally, the NEPTUNE and 
IBM PC-AT should be configured to perform more intelligent types of control 
such as: adaptive, obstacle avoidance, etc. 

Additionally, the following hardware improvements should be made to 
the system. First, it is recommended that the system pressure be increased to 
that specified by the valve manufacturer, and the system should be 
remodeled. It is also recommended that the rack and pinion joint 
transmission gears be redesigned to remove or minimize backlash. Lastly, 
digital shaft encoders should be installed to increase the accuracy of joint 
position feedback. 
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APPENDIX 

OPTIMAL CONTROL PROGRAM 



$include:'labhead.for' 

$include:’master.for' 

^includei'dadio.for' 

dimension chan(6),databuf(6,2500) 

real*4 volt,runt,terml,term2,svltlr,svlt2r,kopt(2,5),deltjl,deltj2 
real*4 stptl,stpt2,endptl,endpt2, lambda, thdotl 
logical flagl,flag2 

c constants 

thlzero = 41 
th2zero = 1300 
thlnint = 2620 
th2nint = 3835 
thlrng = thlnint-thlzero 
th2rng = th2nint-th2zero 
count =1 
runt = 40 
nsweep = runt*20 
dither = 10 

c set channel numbers to be scanned 

chan(l) = l 
chan(2) = 2 
chan(3) = 6 
chan(4) = 7 
chan(5) = 8 
chan(6) = 9 



c reset and initialize labpac 



result = labpac(RESET) 
result = labpac(ATOD,12,AIINIT) 

c initialize dadio board for digital to analog output 

result = labpac(DA,4,AOINIT) 
result =labpac(7, 125, TIMER, SWINIT) 
result = labpac(0,0, AO RAW) 
result = labpac( 1 ,0, AO RAW) 



c read in the starting and ending joint positions 
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400 write (*,*) 'ENTER JOINT 1 STARTING POSITION (DEGREES) 
read (*,*) stptl 

write (*,*) 'ENTER JOINT 1 ENDING POSITION (DEGREES) 
read (*,*) endptl 

write (*,*) TENTER JOINT 2 STARTING POSITION (DEGREES) 
rea( j ^ ^ j stpt2 

write (*,*) ’ENTER JOINT 2 ENDING POSITION (DEGREES) 
read (*,*) endpt2 
write (*,*) 
write (*,*) 

c read in the array of control feedback constants 

write(*, *) ’ENTER THE FILE NAME THAT CONTAINS THE K 
+ MATRIX’ 
write (*,*) 

open (12, file = ”,status = ’old') 
read(12,*) lambda 
do 50 i = 1,2 

50 read (12,*) kopt(i,l),kopt(i,2),kopt(i,3),kopt(i,4),kopt(i,5) 



c move joint 1 and 2 to starting position 

pause 'TURN ON MASTER SERVO POWER SWITCHES’ 

fiagl = . false. 

flag2 = . false. * 

term2 = stptl/90.*thlrng 

iterm2 = ifix(term2) 

100 thl =labpac(l,AIRAW) 
delthl =thl-iterm2 
svltlr = -(de 1th 1*2047. )/thlrng 
svoitl =ifix(svltlr) -fdither 
if(svoltl .It. -2047) then 
svoitl =-2047 
endif 

if (svoitl .gt. 2047) then 
svoitl =2047 
endif 

if(abs(delthl).le. 15) then 
result = labpac(0,0,AORAW) 
fiagl = . true, 
else 

result = labpac(0, svoitl, AORAW) 
fiagl = . false, 
endif 

th2 = labpac(2,AIRAW) 
terml = stpt2*(th2rng/90.) 
iterml =inx( terml) 
delth2 = th2-th2nint-iterml 4- thl 
svlt2r = -(delth2*2047.)/th2rng 
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svolt2 = ifix(svlt2r) + dither 
if (svolt2 .It. -2047) then 
svolt2 = -2047 
endif 

if (svolt2 .gt. 2047) then 
svolt2 = 2047 
endif 

if (abs(delth2) .le. 15) then 
result = labpac( 1 ,0 , AORAW) 
flag2 = .true, 
else 

result = labpac( l,svolt2, AORAW) 
flag2 = .false, 
endif 

if (flagl .and. flag2) then 
goto 200 
endif 

dither = -dither 
goto 100 



c commence optimal control and data acquisition run 



200 write(*,*) 

pause 'READY TO COMMENCE OPTIMAL CONTROL MOVE' 



c convert theta 2 from global to relative 
endpt2 = 180.-endptl + endpt2 



swleft = nsweep- 1 
lstoutl =0. 
lstout2 = 0. 

result = tabpac(l, NSWEEP, 6, CHAN, DATABUF.AISWST) 



c optimal control algorithm with 20th of a second time delay 
c between control sweeps 

300 result = labpac(SWLEFT,AIST AT) 

deltjl =float(databuf(l,count)-thlzero)*90./float(thlrng)-endptl 
deltj2 = float(databufI2,count)-th2zero)*90./float(th2rng) + 90. 

+ -endpt2 

if(((abs(delt]T).le. 1.) .and. (abs(deltj2).le. 1.)) .or. (swleft 
+ .eq. 0)) then 

result =labpac(0,0, AORAW) 
result = labpac( 1 ,0, AO RAW) 
result = labpac(AISWAB) 
stop 
endif 

thdotl = float(databufl( l,count-l)-databuf(l, count))* 1800,/float 
+ (thlrng) 
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svltlr = -kopt(l,l)*deltjl-kopt(l,2)*thdotl-kopt(l,3)*deltj2 
+ kopt(l,4)*lstoutl-kopt(l,5)*lstout2 
svoltl = ifix(svltlr*2047/10. + dither) 
svlt2r = -kopt(2,l)*deltjl-kopt(2,2)*thdotl-kopt(2,3)*deltj2 
+ kopt(2,4)*lstoutl-kopt(2,5)*lstout2 
svolt2 = ifix(svlt2r*2047/10. + dither) 
result = labpac(0, svoltl, AO RAW) 
result = labpac( l,svolt2,AORAW) 
lstoutl =svltlr 
lstout2 = svlt2r 



count = 
swleft = 
dither = 



count + 1 
:swleft-l 
-dither 



go to 300 
end 
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Section IV 



Robot Firefighter Drilling Loads 
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LIST OF VARIABLES 



A 


The cross-sectional area (in 2 ) 


c 


The distance from the neutral axis to the most 
extreme fiber (in) 


T 


Support tube mean diameter (0.9 in) 


E 


Young's Modulus of Elasticity (29,000,000 lbs/in 2 ) 



F r , F r ,F r Resultant forces for an arbitrary point 
x Y on z the support tube cylinder (lbs) 



F l' f 2 f f 3 


Load forces at the penetrator tip (lbs) 


G 


Shear Modulus of Elasticity (lbs/ in) 


I 


The moment of inertia of the cross-section with 
respect to the neutral axis (in 4 ) 


J 


The polar moment of inertia of the cross-section 
relative to its centroid (in 4 ) 


L 21' L 22 


The distance from the centerline of the drill to 
the critical points on the support cylinder surface 
(L 2 i = 2.75 in and L 22 = 6.75 in) 


M 


Bending moment (in-lbs) 


M r ,M r „, 

A v 


Vi Tr7 The resultant moments for an arbitrary point 

on z the cylinder (in-lbs) 


Mi,M 2 ,M 3 


The load moments at the penetrator tip (in-lbs) 


p 


Transmitted force (lbs) 


Q 


The first moment of the area with respect to the 
neutral axis (in 2 ) 


r 


Support tube mean radius (0.45 in) 


r i 


Inside radius of the support cylinder (0.4 in) 


r o 


Outside radius of the support cylinder (0.5 in) 


t 


Radial thickness of the support cylinder (0.1 in) 
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T Torque (in-lbs) 

V Transverse shear force (lbs) 

ec The angle between a force vector and the x axis 

(degrees) 



0 The angle between a force vector and the y axis 

(degrees) 



* Shear strain (in/in) 

* X y,Yy Z ,* zx The shear strain where the first subscript 
denotes the plane face on which the strain acts and 
the second the direction on that face (in/ in) 

^yz ( L 2l) The shear strain at the point B a distance L 2 2 _ 

B from the centerline of the penetrator (in/ in) 



ft The angle between a force vector and the Z axis 
(degrees) 



t. Normal strain (in/in) 



£. e c 

X' y ' z 



The normal strain where the subscript denotes 
the plane face on which the strain acts (in/in) 



A strain measured by a rosette gage element in the 
Y-Z plane (in/in) 

€ "va(L2i) The normal strain at a point A on the cylinder a 
distance L 2 i from the centerline of the penetrator 
(in/in) 



e y, (L 22) 



The normal strain at a point A on the cylinder a 
distance L 22 from the centerline of the penetrator 
(in/ in) 



«y A (L2l) 



The normal strain at a point B on the cylinder a 
distance Lot from the centerline of the penetrator 
(in/in) 



( l 22 ) 

A 



The normal strain at a point B on the cylinder a 
distance L 22 from the centerline of the penetrator 
(in/in) 



^y (^21) 
B 



The normal strain at a point C on the cylinder a 
distance L 21 from the centerline of the penetrator 
(in/in) 
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e yz ( L 2l) The Y_z plane strain at a point B on the 

B cylinder a distance L 21 from the centerline of the 

penetrator (in/ in) 

£ 

z B( L 2l) The norina l strain at a point B on the cylinder a 
distance L 21 from the centerline of the penetrator 
(in/ in) 

Q An angle in the X-Z plane measured from the X axis 

(degrees) 

V Poison's ratio (0.29) 

<r Normal strain (lbs/in) 

, <Ty , d* z The normal stress where the subscript denotes 

the plane face on which the strain acts (lbs/in) 



O’y / - The normal stress acting on the Y plane face at 
A B tne points A, B,and C 



<*^(r, 9 ) The normal stress acting on the Y plane face as a 
function of the variables r and 0 (in/in) 

^ Shear stress (lbs/ in) 









z , *y Xr shear stress at the points 

3 where the first subscript denotes the 
on which the stress acts and the 
direction on that face (in/in) 



A, B, and C 
plane face 
second the 
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I. 



INTRODUCTION 



A. PROBLEM DESCRIPTION 

An aircraft skin penetrator/applicator was developed by 
the USAF as a device for penetrating an aircraft fuselage 
and injecting Halon into the interior of the aircraft to 
extinguish fires [Ref. 1]. The guidelines for the present 
study and other related Naval Postgraduate School work were 
to investigate the requirements for a robot arm to 
manipulate such a tool. The purpose of this work was to 
identify the loads that the support structure would have to 
sustain during a drilling operation with the device. In 
general, the loads are functions of the material character- 
istics' of the drilling surface and the drill parameters, 
i.e., type drill bit, length of the drill, feed rate, and ' 
angle of attack. For this investigation, the drill bit 
parameters were predetermined from previous work at the 
Naval Postgraduate School [Ref. 2 }. Sheet aluminum was used 
to simulate the aircraft fuselage. 

B . APPROACH 

The approach was to experimentally determine the drill 
tip loads by selective support structure strain measure- 
ments. The analytical method used was to resolve the loads 
at the drill tip into six force and moment components about 
an orthonormal coordinate axes system (Figure 1) . These 
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Figure 1. Penetrator Forces and Moments 



83 



loads create a state of strain on the outer surface of the 
support structure which was analyzed by conducting a static 
force/moment equilibrium analysis and a solid mechanics 
analysis. The resulting strain equations were evaluated to 
determine a set of selective strain measurements which would 
facilitate segregation of the individual effects of each 
drill tip load. The analytical expressions for the strain 
at those locations were then used to find analytical 
expressions for the force/moment loads. In this way, the 
measured strains on the support structure were used to 
estimate the six forces and moments at the drill tip which 
are experienced during a drilling operation. 

C . PREVIOUS WORK 

Reference 2 is a previous study of the same problem with 
the assumption that the effects of certain forces and 
moments were negligible (Figure 2) . Notice that only two 



forces 


(Fl 


and 


F 3 ) and one 


moment 


<M 3 ) 


were previously 


assumed 


to 


be 


significant. 


The 


same 


methodology for 


determining 


ana 


lytical expressions 


for 


the force/moment 



components as functions of measured strain quantities was 
used in both studies, but the present study removes this 
restrictive assumption on the load magnitudes. The previous 
work also used an Operational Test and Evaluation drill to 
make measurements, while this investigation used a new drill 
supplied by the NSWC/White Oak laboratory. 
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Figure 2. Penetrator Forces and Moments of Previous Work 
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II. EXPERIMENTAL EQUIPMENT 



A. PENETRATOR ASSEMBLY AND FOUNDATION 

1. Aircraft Skin Penetrator/Applicator 

An AMETEK, Offshore Research and Engineering 
Division model ASP/AA-1 skin penetrator/applicator was 
utilized (Figure 3) . For this work the carrier and bottle 
assemblies and the Halon ball valve were removed to 
facilitate installation on a milling machine bed. The. drill 
was designed to be powered by compressed air carried by the 
operator. In the present work, a source of compressed shop 
air (165 psig) was used to power the drill. 

2 . Instrumented Support Fixture 

A support fixture was manufactured from mild steel 
which consisted of a central cylinder with top' and bottom 
support plates welded to its external radius. The 

penetrator was attached to the top plate with eight 
self-locking helicoil fasteners. The bottom plate was 

bolted to a milling machine bed. Five strain gages were 

mounted on the cylinder surface to facilitate the 
measurement of the cylinder deformation during a drilling 
operation. 

3 . Milling Machine 

A Milwaukee model H milling machine, with mill 
removed, provided a movable foundation for the penetrator. 
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Figure 3. Experimental Equipment Left View 





Two Superior Electronic synchronous/stepping motors, 72 RPM, 
120 volts and 60 hertz, were installed on the milling 
machine's transversing gear assemblies and were controlled 
by 3-position switches. These motors provided for movement 
along the Y and Z directions at a constant speed of 14 
inches/minute . 

B. SIMULATED AIRCRAFT FUSELAGE 

1 . Test Panels 

As in Reference 2, the aircraft fuselage was 
simulated by test panels constructed of 7075-T6 aluminum 
0.050 inch sheets. The panels were rolled to the shape of 
the support frame. 

2 . Support Frame 

The support frame with a test panel installed is 
shown in Figure 4. The frame was bolted to the floor and 
leveled. Due to the manner in which the support frame was 
constructed, the minimum angle of attack achievable, without 
penetrating the support frame, was 8.48 degrees between the 
drill axis and the sheet aluminum. 

C. STRAIN MEASUREMENT/RECORDING EQUIPMENT 

1 . Strain Gages 

Linear gages were installed at critical points on 
the support structure. These gages were Micro-Measurement 
Division series CEA-06-125UN-350 constantan strain gages 
with a gage factor of 2.065 +/- 0.5% , and a resistance of 
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Figure 4. Experimental Equipment Right 



3 50 ohms +/- 0.3%. A rosettee strain gage was used to 
determine shear force; it was a Micro-Measurement Division 3 
element 45 degree rectangular stacked rosette, series 
WK-06-120WR-350, superimposed K-alloy, strain gage with a 
gage factor of 2.065 and a resistance of 350.0 ohms +/- 
0.4%. 

2 . Recording Equipment 

Recording of the strain gage outputs was 
accomplished with an Astro-Med, Inc. 8 channel strip chart 
recorder via an Ectron amplifier assembly. 
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III. ANALYTICAL METHOD 



A. BACKGROUND 

1. System Model 

The system model was a rigid body assembly which 
consisted of a pneumatic drill mounted on a rigid, movable 
support structure. The support structure was a thin-walled 
tube (diameter >> thickness) welded to support plates top 
and bottom and fixed to a rigid milling machine bed. 

2 . • Theoretical Approach 

The first step was to conduct a static analysis of 
the system to determine the force/moment resultants at any 
arbitrary point on the support tube as a function of the 
drill tip loading. Following this, a solid mechanic 
analysis was conducted to determine the particular effect of 
each force/moment component. Using the theory of 

superposition, the six effects were added to determine the 
net effect of the force/moment components. These were then 
converted to functions of drill tip loading by utilizing 
the equilibrium equations of the static analysis. The next 
step was to apply the fundamental triaxial stress-strain 
relationships. This provided analytical expressions for the 
strain at any arbitrary point on the tube as a function of 
the loads at the drill tip. The final step was to evaluate 
the reverse case; that is, to determine analytical 
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expressions of the force/moment loads of the drill tip as 
functions of the strains at an arbitrary point on the tube. 
3 . Theory of Superposition 

The theory of superposition for linear systems 
states that the net stresses of a system loaded with a 
combination of forces and moments can be analyzed by 
determining the resulting stresses of each load and adding 
the like stress components. Therefore, an analysis of the 
normal strains and shear strains resulting from the reaction 
loads of Figure 1 at an arbitrary cut on the tube was 
analytically evaluated for each load and then combined to 
determine the net effect. 

B. STATIC SOLID MECHANICS ANALYSIS 
1 . Static Analysis 

A free body diagram of the system with static 
loading is presented in Figure 5. The static equations are: 



I F x “ 0 : 

I F y = 0: 

X F z - 0: 
£ M r = 0 : 



F 1 + Fy- — 0 Fy- — “F-J 



F 2 + F r - o -> F r - -F 2 



F z + F r y “ 0 F r z ~ ~ f 3 

A A AAA A 

(L 2 j - L^J&) X (F X i+F 2 j + F 3 k) + (M!+M rx )i 



+ (M 2 +M r )j + (M 3 +M r )k = 0 

y z 



Equating the like components equal to zero yields: 
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Figure 5. Penetrator Free Body Diagram 
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i components: L 2 F 3 + + M r = 0 

A 



+ M ry “ “ ( L 2 F 3 +L 1 F 2 +M l) 



j components: “L 1 F 1 + M 2 + M r 



0 - % = L 1 F 1 - 



k components: -L 2 F 1 + M 3 + M r 



z 



M r z “ l 2 f 1 “ m 3 



2 . Load-Stress Fundamental Relationships FRef. 31 

For axial tension, the normal stress is given by: 



P = the transmitted force 
A = the cross-sectional area (A = irtd) 

For axial bending, the normal stress is given by: 



M = the bending moment 

c = the distance from the neutral axis to the 
most extreme fiber (c = r Q ) 

I = the moment of inertia of the cross-section 
with respect to the neutral axis 



For transverse loading, the shearing stress is given 



<r 



p 

A 



where: 




where : 



by: 
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where 



r 



VQ 

Jt 



V = the transverse shear force 

Q = the first moment of the area of the semi- 
circular section with respect to the 
neutral axis 

rn _ 2, 3 3,, 

[Q - T (r -r. )] 

3 o l 

t = the radial thickness of the cross-section 

J = the polar moment of inertia of the cross- 
section relative to its centroid 

<j = ftd 3 ) 

For torsion of thin-walled tubes, the shearing 
stress is given by: 






Tr 

J 



where: 

T = the applied torque 
r = mean radius 

J = the polar moment of inertia of the cross- 
section relative to its centroid 

<j = i-td 3 ) 
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3 . Stress Analysis 

The state of stress at a point on the tube is 
described by nine values of stress. The following is an 
analysis of each reaction load to determine their 
contribution to these stress values; 

For M r , Figure 6 depicts the loading condition, 
z 

with the resulting stress distribution. Notice that the 
stress at A and C is zero, and at B is maximum. These 
critical points will be important in all the stress and 
strain relationships which follow. The stress distribution 
is given by: 



<sMr,d) 



M r cos 0 
r z ° 



For M 



r„ ' 



Figure 7 depicts the loading condition, 



with the resulting stress distribution. The 
distribution is given by: 



stress 



«* y Cr,0) 



M r sin -0 
r o 
x 



For M 



r , Figure 8 depicts the loading condition, 

The stress 



with the resulting stress distribution, 
distribution is given by: 



X. = 



M r 
r 

y . _ 
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Figure 6. 






* y 




M r Load and Stress Distribution 

z 
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Figure 7. M r Load and Stress Distribution 
x 
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X 



Point C 



~A 






Ti 



Figure 8. 




Load and Stress Distribution 
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For F r , Figure 9 depicts the loading condition, 
z 

with the resulting stress distribution. The stress 

distribution is given by: 



r * 




Q 



Jt 



For F r , Figure 10 depicts the loading condition, 
with the resulting stress distribution. The stress 
distribution is given by: 




For F 



r , Figure 11 depicts the loading condition, 



with the resulting stress distribution, 
distribution is given by: 



The stress 



F 

r 




Now, utilizing the theory of superposition the net 
stress at points A, B, and C can be determined. 
Additionally, utilizing the static equilibrium equations of 
Section II.B.l, the following six stress equations are 
produced: 



Point A: <r 

^A 





- (L 2 F 3 +L 1 F 2 +M 1 ) . F 2 
I A 




F rx Q _ - (L l F r M 2 )r , F 1 Q 
Jt J Jt 
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Point C 
















Figure 9. F r Load and Stress Distribution 
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Point B 
^V~7 



Point A 








F 

rx 



Point C 





i i 




M 

-k. 



Figure 10. F r Load and Stress Distribution 
x 
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Figure 11. F r Load and Stress Distribution 

y 
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Point B 



r B 



yz 



B 



M r 
r z ° 



-M r 
r 

y 



F 

r 

_Z 

A 



Jt 



(L 2 F rt )r o F 2 

I + A 



-(L^-M^r F 3 Q 
J + ~Jt~ 



Point C: 



-M r F 
r o r 

£ Z = 

I A 



(L 2 F 3 +L 1 F 2-«I ),: o , F 2 
I A 



yx-c 



M r 
r 

_y 

j 



jt 



(L 1 F 1 -M 2 )r F X Q 
J + “jt 



4 . Stress-strain Fundamental Relationships 



e 


= ho - 


v(o 


+ a )] 


X 


E x 


y 


z 


£ 


= ho - 


v (a 


+ a ) ] 


y 


E y 


x 


z 


£, 


= ito - 


v(a 


+ a )] 


-7 


E 2 


X 


y 


*xy 


li 

< 

O 








n 

o 






^zx 


= T zx/ G 







5 . Determination of Strains 

As stated in the theoretical approach, Section 
II. A. 2, the next level of analysis was to determine 
analytical expressions for the strains at points A, B, and C 
as functions of the loads. This was accomplished by 
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substituting the expressions for stresses of Section II. B. 3 
into the fundamental triaxial stress-strain relationships of 
Section II. B. 4. Table I summarizes the results. Appendix A 
contains the Fortran computer code for these equations. 

6 . Determination of the Drill Tip Loads 

Figure 12 shows the locations of the strain sites of 
interest. Point A is a point on the outer surface of the 
tube in the Y-Z plane. L 2 i is the distance from the 
centerline of the drill to any arbitrary Z-X plane cut of 
the tube (e.g. , horizontal cut) through point A. The normal 
strain at point A is a function of L 2 i only and is given by: 




( L 2i ) 



l f -(L 2i F 3 +L l F 2 +M l )r o F 2, 
E l I A J 



Therefore, taking the difference of two normal strains whose 
Y coordinates are L 2 i = L 2 2 and L 2 i = L 2 i, such that L 22 > 
L 2 1 , yields 



e (L„„) 

y A 22 




(L 21> 



( L 22 +L 21 )r o F 3 
El 



Solving for F 2 yields 



F 3 = 



[C y A (L 2 2) ~ e y A (L 21 )]EI 
( " L 22 +L 21 )r o 



( 1 ) 
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STRAIN EQUATIONS 
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a. Side View (Y-Z plane) 




b. Front View (X-Y plane) 



Figure 12. Strain Measurement Locations 
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Alternatively, using the same rationale for point B 
and taking the difference of two normal strains with 
coordinates Y = L 2 2 and Y = L 21' yields 



e y B (L 22 ) £ y F (L 21 ) 



f B 



(L 22' L 21 ) r o F 1 
El 



Solving for yields 



t«y < L 22> -V <L 2 i’ 1EI 

B 



F i " 



B 



(L 22 " L 21 ) r o 



( 2 ) 



At points A and C, the M r bending moment normal 
stress components are equal but in opposite directions. 
Therefore adding the normal strains at these points yields 



<y A (L 21> + e y c !L 22 1 



2F 2 

EA 



Solving for F 2 yields 



p = 
‘2 



[ V (L 21> +e y c ,L 21 !1EA 



( 3 ) 



With the three load forces known, the moment loads 
can be determined directly from the strain equations. 
Solving the normal strain at point A at L 2 i for yields 



.. _ Irr, It n , F 2 (L 21 F 3 +L l F 2 )r o, 

\ ~ r _I " Ee y. (L 21 } A I ] 

o 



(4) 
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Solving the normal strain at point B at L 21 for M 3 yields 



F 9 L F. r 

Ms = ^[-Ee (L 21 , + -l + JiA°] 

O B 



(5) 



Solving the shear stress at point B at L 2 i for M 2 yields 



F 3 Q Lj^r 



*2 * = [Gy y ZB (L 21> - -5t + — 1 



( 6 ) 



7 . Previous Load Determinations 

As stated earlier, the previous work assumed the 
principal loads were F]_ , F 3 , and M 3 and the other loads 
were negligible. The further assumption that F^ did not 
contribute to the bending moment about the Z axis was also 
made. Thus the strain equations of Table I reduce to: 



-r 



% (L ) 



J 21 ‘ 3 ‘o 
El 



-M 7 r 

V L 2i } = -zr 



1 " L l F l r F 3 Q 
*y (L 21 } = + 

Z B 



Solving for the loads yields 



F, = f - ~ E * € (L_. ) 

3 L 21 r o ^fl 21 



- — 






(7) 



( 8 ) 
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( 9 ) 



J F 3 Q 

Vi = - § [g v 2 i> ■ ' 

In a later part of this report, the results of equations 
(1) , (5), (6) and equations (7), (8), (9) will be compared 

in order to evaluate their predictions. 

C. STRAIN GAGE PLACEMENT 

1 . Strain Gage Type and Location 

In order to solve equations (1), (2), and (3) the 

axial strain at five locations is required. Additionally, 
the shear strain is required at one location for equation 
(6) . Table II summarizes the strain gage types and 
locations, and Figure 12 shows their locations. 

2 . Determination of the Shear Stress 

Although the normal strains can be measured directly 
from the gages, additional computations are required to 
determine the shear stress. 

Figure 13 depicts a 45 degree rectangular rosette 
strain gage. The applicable strain equation is [Ref. 3]: 





Substituting this expression into the expression for M 2 
yields 



t F 3 Q L l F l r 

M2 = = tG[2E yz B (L 21 ) " £ z b (L 21 ) _£ y B (L 21 )1 " "jt + ~ J~ ] 
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TABLE II 



STRAIN GAGE TYPES AND LOCATIONS 



Point 


Moment 

Arm 


Type 


Quantity 

Measured 


A1 


l 2 1 


Linear 




A2 


l 22 


Linear 


Cy 


B1 


l 21 


45 degree 

rectangular 

rosette 




B2 


L 22 


Linear 


€ y 


Cl 


L 2 1 


Linear 


Cv 




Figure 13. 45 Degree Rectangular Rosette 



Table III summarizes the force/moments equations for this 
work and Table IV for the previous study. Appendix B 
contains the Fortran computer code for these equations. 

D. CALIBRATION PROCEDURES 
1. Approach 

Static loads were applied to the support structure 
to evaluate whether the analytical expressions of Table I 
correctly predict the measured strains at locations A, B, 
and C. Figures 14, 15, 16, and 17 depict the methods by 

which the static loads were applied to the instrumented 
drill support structure. The static loads consisted of 
round stock cut into approximately ten pound segments. Each 
segment weight was measured using an Instron Corporation 
Material Testing System 100 pound test cell. For this test 
cell the accuracy of the measurement was +/- 1.0 % of the 
indicated load, or 0.02 % of the test cell capacity (100 

lbs) . 

The weight holding apparatus depicted in the static 
loading figures was of sufficient strength to hold the 
required testing loads. However, due to misalignment in the 
support structure, the loads were not applied exactly . along 
the desired coordinate directions associated with the 
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TABLE III 



LOAD EQUATIONS 



F, = 



|e y B (L 22> - e y B <L 21 )IEI 
<L 22 " L 21 )r o 



F- = 



l£ y A (L 21> - e y c (L 21 )1EA 



F, = 



[e y < L 22> - V (L 21 )1EI 

/i H 

(~ L 22 +L 21 )r o 



_ F_ (L_ n F_ +L 1 F_)r 

“l - ~ L 21 + - A I 1 

O -'A 



f 3 q 



f-U = — (G[2e (L„,)-c - ~^r + 

2 — yz„ 21 z„ 21 v„ 21 Jt 



B -B 



*2 ■ r‘- E£ y < L 21> + 4 + -^T^' 

O -B 



TABLE IV 

LOAD EQUATIONS OF PREVIOUS WORK 

F = -^L. e ' ■ . 

3 r o L 21 ^ 

“3 = -¥ M 3 
o 

Vl “ - § (G f 2£ yz B (L 21 ) - e z B (L 21> - £ y B (L 21 )J 



* 




Figure 14. F 2 Static Loading 







Figure 15. Static Loading 
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Figure 16 . F3 Static Loading 
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#•%&# 






support structure. Figure 18 illustrates the problem 



encountered when a load was attempted along the Z-axis. To 
overcome this obstacle, an angle position indicator was used 
to measure the angles between the orthonormal axes and the 
line of action of the static loads. Only two of the three 
angles were measurable and the following relationship was 
used to determine the third angle: 



The component of the static load acting on each of the 
orthonormal axes was then determined using the following 
relationships : 



As can be seen from Table III, each of the forces is 
equal to a strain quantity times a constant, i.e. 



1 




F]_ = F cos o( 



= F cos p 



F 3 = F cos & 



F 1 = Uy B ( L 22) - £ y B ( L 2l) ] Constant 



F 2 = [ £ y A (L 2 i) + e y c (L 2 i) ] Constant 



f 3 = C € y A ( L 22) " ^^ 21 ) ]*constant 
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Y 




Figure 18. Static Force Line of Action 

Therefore, in comparing the measured versus calculated 
values, these strain quantities were used to determine the 
accuracy of the analytical strain equations. 

2 . Static Loading Results 

Figure 19 shows the plot of [ £y (L 2 2 ) “ e y B ( L 2l) ] • 
both measured and calculated, versus an applied F^ static 
load of 30.74 to 99.65 pounds. As can be seen from the 
figure, there is a negligible difference between the 
measured and calculated values. Figures 20 and 21 show the 
plots of the F 2 and F 3 component forces with their 
respective strain quantities. The substantial difference 
shown in the later figures will be discussed below. 
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EYR22-LYH21 (MICRO-STRAINS) 

0.0 30.0 60 0 WOO 120.0 150.0 180.0 210.0 240.0 



FI STATIC TEST 




Figure 19. Load, F^ Plot 
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EVA:.'! + r.YC21 (MICRO-STRAINS) 

(i O 8 0 10.0 12.0 14.0 16.0 



FI STATIC TEST 




Figure 20. Load, F 2 Plot 
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EYA22-EYA21 (MICRO-STRAINS) 

e 0 12.0 15.0 18.0 21.0 24.0 



FI STATIC TEST 



LEGEND 

a CALCULATED. 
"5 MEASURE IT 



i 



© 

d ' 



© 
CO * 




z 0 



o 

d _L 



0.0 1.0 2.0 3.0 

-F3 (LBS) 



4.0 



5.0 



Figure 21. Load, F 3 Plot 
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Figure 22 shows the measured versus calculated 
strain plot for an F 3 load of from 30.77 to 198.66 pounds. 
Again the measured and calculated values of the principal 
strain quantity correlate well. Figures 23 and 24, for the 
transverse force components, exhibit the same disporportion- 
ate errors as was found with the F 3 static load. 

Figures 25 and 2 6 show the plots for major X and Z 
axis loading. This loading was accomplished by placing the 
support fixture at an angle of 44.5 degrees with respect to 
the load line of action. Therefore, a static load of 30.77 
to 198.66 pounds produced F 3 and F 3 loads of 22.83 to 147.42 
pounds and 20.61 to 133.06 pounds respectively. Again the 
strain equations produced results very close to the measured 
values with the exception of the transverse F 2 force, as 
shown in Figure 27. 

The results of applying an F 2 static load are shown 
in Figure 28. F 2 loading was accomplished by placing 
weights of 10.0 to 100.02 pounds directly on top of the 
support fixture. In this test, there were considerable 
differences in the measured and calculated values. 

3 . Sources of Error 

The errors encountered with the transverse force 
components and F 2 loading are attributed to the nonsymmet- 
ric shape of the tube, the nonrigidity of the tube, and the 
nonhomogeneous tube material. ' These factors all contribute 
to change the centroidal axis location. 
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— (EYA22- EY2J) (MICRO-STRAINS) 

0 180 0 210.0 300.0 360.0 420.0 480.0 



F3 STATIC TEST 




Figure 22. F 3 Load, F 3 Plot 
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n (MICRO-STRAINS) 

6.0 <10 12.0 > 5.0 18 0 21 0 24.0 



F3 STATIC TEST 




Figure 23. F 3 Load, F^ Plot 
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rVA2l + EYC21 (MICRO-STRAINS) 

■ 1.0 6.0 8.0 10.0 12.0 



F3 STATIC TEST 




Figure 24. F 3 Load, F 3 Plot 
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KYUi^-EYDSl (MICRO-STRAINS) 

OO 45.0 00 0 135 0 180.0 225.0 270.0 315.0 360.0 






FI AND F3 STATIC TEST 




0.0 30.0 



60.0 90.0 



120.0 150 



FI (LBS) 



Figure 25. 



Fi and F 3 Loads, F^ Plot 
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-(EYA22-EYA21) (MICRO-STRAINS) 

0.0 40.0 80.0 120.0 160.0 200.0 240.0 280.0 320.0 






FI AND F3 STATIC TEST 




Figure 26. F 3 and F 3 Loads, F 3 Plot 



128 



EYA2H KYC21 (MICRO-STRAINS) 
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FI AND F3 STATIC TEST 
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Figure 27. and F 3 Loads, F 2 Plot 
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-(tVALM+EU'21) (MICRO-STRAINS) 

D 12 0 16.0 20 0 24.0 28.0 



F2 STATIC TEST 



o 




Figure 28. F 2 Load, F 2 Plot 
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As can be seen in Figures 6 and 7 the location of 
the netitral axis for the bending moments coincides with the 
centroidal axis. Also, the magnitude of the stress is 
directly proportional to the distance from the neutral axis, 
and the stress distribution is symmetrical about it. These 
facts operate with the non-ideal, real tube discussed above 
to create a changing location of the centroidal axis and, 
consequently distort the stress distribution. As shown in 
the transverse force component plots, the change in 
centroidal axis has a substantial impact on the stress 
distribution and, consequently, the reliability of the 
strain equations to determine the loading. But due to the 
relative magnitudes shown in the preceding figures, these 
errors were considered negligible in predicting the dynamic 
loads during a drilling operation. 

Similarly, the stress distribution resulting from an 
F 2 loading, Figure 11, was nonuniform with a change in the 
centroidal axis location, yet agreement was good. 
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IV. DRILLING OPERATION AND RESULTS 



* 



A. EXPERIMENTAL PROCEDURES 

The strain gage amplifier assembly and strip chart 
recorder were energized and given a one hour warm-up period 
prior to commencement of recording data. The air compressor 
was started one-half-hour prior to the drilling operations 
to allow sufficient time to pressurize the accumlator to 167 
psig. After the one hour warm-up period, the amplifiers 
were calibrated to 4.84 volts per 10,000 micro-strains. The 
penetrator trigger was taped in the full open position and 
air flow was controlled by an in-line valve. The penetrator 
was positioned via the milling machine's two motors to 
permit the penetrator an angle of attack of 8.48 degrees and 
a three diameter separation between drill hole centers. 

At the commencement of a drilling operation, the in-line 
air valve was fully opened and the Z-axis motor of the 
milling machine energized to a bed speed of 14 inches per 
minute. All systems were stopped after the penetrator' s 
drill had completed drilling the hole. 

The strip . chart recordings were then converted to 
micro-strains. A sample of the strip chart readings is 
shown in Figure 29. Since the conclusions of this work are 
to be used in the design of the penetrator support 
structure, a worst-case set of data was needed. Therefore, 
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Figure 29. Sample of Strip Chart Recording 



a french curve was used to curve-fit the maximum values of 
the recorded strains. To accurately reflect the continous 
time effects of the drilling, data at every 0.2 seconds was 
used. The micro-strain data was then entered into the 
static analysis Fortran program of Appendix B, which is the 
computer code equivalent of the force/moment equations of 
Table III and Table IV. 

B. AIR SUPPLY PRESSURE 

To obtain viable data, it was essential that the 
penetrator successfully complete the hole drilling without 
binding. The previous work used an operating pressure for 
the penetrator of 100 psig [Ref. 2]. In the present work, 
a pressure regulator was initially installed in the air line 
to reduce the accumulator pressure from 167 to 100 psig. 
With 100 psig supply pressure, the penetrator was 
unsuccessful in completing a drilling operation without 
binding. The supply pressure was increased in steps of 25 
psig from 100 to 150 psig and no run was achieved without 
experiencing binding. Therefore the pressure regulator was 
removed from the system and runs were achieved without 
binding with pressures of 158 arid 164 psig. 

As drillings were accomplished, the accumulator supply 
pressure decreased from 167 to 150 psig, at which point the 
compressor would recharge -the accumalators . During the ex- 
perimental process it was noted that any attempts to drill 
with supply pressures of less than 158 psig binding occured. 
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C. DRILLING DATA RESULTS 

Three drilling operations were conducted with supply 
pressures of 158, 164, and 164 psig. For the plots, the 

drillings are differentiated by their duration times as 
illustrated in Table V. 

TABLE V 

LOAD RESULTS DESIGNATION 



Run 

No. 


Pressure 

(PSIG} 


Duration 
Time (sec} 


Functions 


1 


158 


6.0 


Fx<6.0) ,F 2 (6.0) ,F 3 (6.0) , 
Mi (6.0) ,M 2 (6.0) ,M 3 (6.0) 


2 


164 


7.2 


F X (7.2) ,F 2 (7.2) ,F 3 (7.2) 
Mi (7.2) ,M 2 (7.2) ,M 3 (7.2) 


3 


164 


10 . 0 


F 1 (10.0) ,F 2 (10.0) ,F 3 ( 10.0) , 
M X (7.2) ,M 2 (7.2) ,M 3 (7.2) 




1. Forces 








The force loads consist of a thrust component, F 3 , 


acting along the axis of 


the drill and two transverse 


components, F^ 


and F 2 , acting in a plane perpendicular to 


the 


axis of the 


drill. 





As can be seen on Figure 30, the force F 3 was always 
in the positive direction, i.e., pushing against the drill. 
The graphs indicate the values of F 3 are fairly consistent 
no matter what the drilling time was. For instance, the 
maximum values were 89.34, 88.47, and 98.01 pounds for the 
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three runs. Also, the data shows that F 3 increases to a 
peak as the drill penetrates, then drops in magnitude prior 
to increasing to a final maximum peak. The data relates to 
the geometry of the drill tip (Figure 31) as the drill 
penetrates. The first peak is caused by the knuckle, the 
second by the trailing edge. 

The transverse forces on the other hand, vary in 
direction from one run to the next (Figures 32, 33). The 
direction change is attributed to the movement of the drill 
tip prior to stablizing in a drilling configuration; i.e., 
seating into the test panel. If the initial biting movement 
of the drill tip is referred to as "walking," then the 
direction and distance of the walking contribute to both the 
magnitude and direction of the transverse loads. Once a 
stable configuration was achieved in a run the force 
direction for that individual run did not change. In run 1, 
the direction of the drill was down and to the left on the 
surface of the test panel. While for run 2, the direction 
was up and to the left; for run 3, the direction was down 
and to the right. The maximum values were for F^ : -24.28, 
-13.66, and 44.67 pounds; and for F 2 : -50.76, 42.35, and 
-42.31 pounds. The maximum values for F^ occured at the end 
of the drill run and for F 2 at the beginning of the run. 
Therefore, F 2 was the governing transverse force in 
determining the stablized configuration. The F^ force was 
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Figure 30. 
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Figure 31. Drill Bit Schematic 
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FI(LUS) 



Fl: X AXIS FORCE 



O 




Figure 32. Drilling Results for F 1 Plot 
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F2(I.BS) 



F2: Y AXIS FORCE 




Figure 33. Drilling Results for F 2 Plot 
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characterized by increasing magnitude as the drill’s 
trailing edge was approached. The F 2 force exhibits sharp, 
momentary increases in magnitude followed by relatively long 
periods at constant value. 

2 . Moments 

The direction and magnitude of the moments were 
coupled to the effects of the transverse load seating 
conditions. For run 1, the M^ moment about the X axis was 
positive with a maximum value of 1137.88 inch-pounds (Figure 
34), and for run 2, M^ was negative with a maximum value of 
-996.75 inch-pounds. While for run 3, the value of 
varied from -37.53 to 1003.91 inch pounds. Each time, the 
maximum values occured at the beginning of the drilling. 

The M 2 moment about the Y axis oscillated between 
negative and positive values in all three runs (Figure 35) . 
The ranges of values were as follows: 

Run 1: -427.45 to 72.28 inch-pounds 
Run 2: -147.11 to 47.99 inch-pounds 
Run 3: -309.63 to 803.71 inch-pounds 

The peak values. were at the trailing edge of. the drill bit. 

The M 3 moments about the Z axis were also 
oscillatory, but were, for the most part, positive in nature 
(Figure 3 6) . This positive nature correlates to the fact 
that binding did not occur for these runs. The range of 
values obtained were: 
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Ml: X AXIS MOMENT 
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Figure 34. Drilling Results for Plot . 
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M2(1N-LBS) 

- 450.0 - 290.0 - 130.0 30 0 190.0 350.0 510.0 C 70.0 830.0 



M2: Y AXIS MOMENT 
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Figure 35. Drilling Results for M 2 Plot 
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M3: Z AXIS MOMENT 
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Figure 36. Drilling Results for M 3 Plot 
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Run 1: -3.37 to 63.23 inch-pounds 

Run 2: 0.0 to 67.55 inch-pounds 

Run 3: -18.03 to 160.52 inch-pounds 

As with M 2 , the peak values were at the drill's trailing 
edge. 

D. COMPARISON WITH PREVIOUS WORK 

Figures 37 through 45 are plots which compare the 
results of this work to the forces/moments computed using 
the equations of the previous work in Table IV. In the 
plots, various quantities are differentiated by a PW 
preceeding values of the previous work. The values of F 3 , 
Li*Fi, and M 3 correlated fairly well for the two studies and 
thus complement each other's validity. 

The differences in the F 3 force are due mainly to the 
assumption in the previous work omitting the M 3 and L 3 *F 2 
bending moments which contribute to the state of strain at 
point A. 

Similarly, the twisting moment at point B is not only a 
function L]_*F^ but also the moment M 2 which, depending on 
its direction, adds or subtracts from L^*F 1 . The normal 
strain at point B, is not only a function of M 3 , but also 
the bending moment caused by L 2 *F 1 . 

E. ERROR ANALYSIS 

For the strip charts, a one-half of a division error 
equates to the following accuracies for the forces: 
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Figure 37. Current Versus Previous Work: F 3 Plot 1 
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F3(LBS) 
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F3: Z AXIS FORCE (M1=M2=F2=0) 
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Figure 38. Current Versus Previous Work: F 3 Plot 2 
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F3: Z AXIS FORCE (M1=M2=F2=0) 




Figure 39. Current Versus Previous Work: F 3 Plot 3 
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L1*F1 (IN-LBS) 

- 450.0 - 290.0 - 130.0 30.0 190.0 350.0 510.0 670.0 630.0 



L1*F1:Y AXIS MOMENT(M1=M2=F2=0) 
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Current Versus Previous Work: 



M 2 Plot 1 
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Figure 41. Current Versus Previous Work: M 2 Plot 2 
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Figure 42. Current Versus Previous Work: M 2 Plot 3 
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M3: Z AXIS MOMENT(M1=M2=F2=0) 
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Figure 44. Current Versus Previous Work: M 3 Plot 2 
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Figure 43. Current Versus Previous Work: M 3 Plot 1 
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M 3(1 X— LBS) 
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M3: Z AXIS MOMENT(M1=M2=F2=0) 




Figure 45. Current Versus Previous Work: H 3 Plot 3 
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FI: +/- 2.17 pounds 

F2: +/- 42.31 pounds 

F3 : +/- 8.67 pounds 

The results of the static tests of the equations yielded 
the following maximum error percentages for measured versus 
calculated values: 

FI: +/" 0.4 % 

F2 : +/- 10.2 % 

F3 : +/- 1.1 % 
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V. CONCLUSIONS 



An investigation of the requirements for a robot arm to 
manipulate an aircraft skin penetrator/ applicator was 
conducted. The penetrator was a USAF-developed device which 
was designed to penetrate the fuselage of an aircraft and 
inject a gaseous fire fighting agent (Halon) into the 
interior of the aircraft to extinguish fires. Experiments 
were designed to identify worst-case drilling loads, and a 
test apparatus was used to mimic the drilling. Test data 
were collected on the effects of drilling through sheet 
aluminum alloy. Forces and moments were measured by strain 
gages during actual drilling. 

The test data indicates that the loads at the drill tip 
are significantly influenced by drill walking prior to 
stabilizing and penetrating of the drilling surface. This 
is attributed to significant transverse load forces (F^ and 
F 2 ) which act in a plane perpendicular to the drill axis. 
These transverse forces act with the length of the 
penetrator to exert significant moments on the support 
structure. Significant factors in the walking phenomenon 
include: the age and condition of the drill tip; the type 
of material to be drilled; and the angle of attack between 
the drill axis and the drilling surface, all of which 
contribute to the direction and the magnitude of the drill 
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walking prior to penetration. The drill walking in turn 
determines the directions and magnitudes of the transient 
loads . 

The following is a summary of the maximum forces and 
moments for each of the runs recorded: 





Run 1 


Run 2 


Run 3 


F 1 (lbs) 


-24.28 


-13.66 


44 . 67 


F 2 (lbs) 


-50.76 


42.35 


-42.31 


F 3 (lbs) 


89.34 


88.47 


98.01 


Mi (in-lbs) 


1137.88 


-996.75 


1003.91 


M 2 (in-lbs) 


-427.45 


-147.11 


803.71 


M 3 (in-lbs) 


63.23 


67.55 


160.52 


From this summary 


, a negative 


F 2 results 


in a positive M^ 



moment on the drill tip where the reverse case is true for 
positive F 2 . As for , its direction effects the resulting 
moment M 2 . The direction of M 2 is positive when F]_ is 
positive and vice versa. With regard to M 3 , F^ is the 
dominant transverse force, but not to the extent it was for 
M 2 , this is true because the moment arm in the M 3 case is 
the distance from the centerline of the drill to a point on 
the support structure which is much less than the length of 
the penetrator in the M 2 case. But as can be seen above, 
the magnitude of M 3 is significantly higher for run 3 where 
a positive F^ was realized vice the other runs when a 
negative F^ resulted from the drill walking. 
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With regard to the previous study of Reference 2, the 
effects of drill wear are an increased axial force and/or 
torque requirement. The location and nature of the wear 
governs which symptom dominates. Chisel edge wear results 
in increased axial force, land wear results in increased 
torque, and wear at the corners results in increases of both 
[Ref. 4]. The thrust (F 3 ) and torque (M 3 ) of the previous 
work for an 8 degree angle of attack were of the magnitudes, 
respectively, of 280 pounds and 475 inch-pounds for 100 psig 
supply air, and 205 pounds and 255 inch-pounds for 140 psig 
supply. Therefore, the combination of the different wear 
effects discussed above combined with an increase of supply 
air pressure to 160 psig, would explain the higher thrust 
and torque required for the worn drill used in the previous 
study and the lower requirements determined in this study. 
This was also visually evident in the resulting sharp, clean 
holes for the current study and the rough, jagged holes 
drilled with the worn drill. Another contributing factor to 
the differences in the torque moment, is the fact that the 
previous study did not deduct the bending effects caused by 
F 3 times the length from the centerline of the drill to the 
strain gage location. Similarly, the previous assumption 
that there was no M 2 at the drill tip, just twisting of the 
support structure caused by F 3 times the length of the 
penetrator, disallowed any correlation between the two 
studies with respect to F 3 and M 2 . 
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VI. RECOMMENDATIONS FOR FUTURE WORK 



Further investigation for two specific areas are 
required. The first is the phenomenon of the drill walking 
and the resulting effects it has on the drill stabilizing 
and penetrating. The drill has a possible choice of four 
quadrants in which to move. These would result from the 
following combinations of transverse forces: 

Fi and F2 
F^ and -F2 
-F]_ and F2 
-F]_ and — F2 

Drilling data should be obtained for all the combinations 
and analyzed for drill support effects and repeatability. 

The second area of concern is the life cycle of the 
penetrator drill. The thrust and torque requirements are 
both a function of drill wear. In view of the fire fighting 
mission of the device, a definitive useful life cycle for 
the drill bit should be identified. - - . 

With regard to the experimental procedures, the 
methodology of this study is sufficient but improved data 
recording is required for further work. A process for 
directly converting the continuous time output of the strain 
gages into an input for the strain-to-load conversion 
computer program should be developed. 
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APPENDIX A 



STRAIN EQUATIONS FORTRAN COMPUTER CODE 

XXX XXX X XXX X XXXX X X X XX X XXXXXXX XX XXX X XXXXX XXXXXXX XXXX X XXXX XXXXXXXXXXXX X X • 

* LAWRENCE, D . A . 

* DRILL1 PROGRAM CLACULATES THE RESULTING STRAINS GIVEN THE APPLIED 

* MOMENTS AND FORCES AT THE DRILL TIP. 
x VARIABLES USED ARE: 

X TH: RADIAL THICKNESS (IN) 

x RO: OUTSIDE CYLINDER DIAMETER (IN) 

x RI: INSIDE CYLINDER DIAMETER (IN) 

X E: YOUNGS MODULUS OF ELASTICITY (LBS/INXX2) 

x P 0 1 : POISON RATIO 

x PI 

x G: SHEAR MODULUS OF ELASTICITY (LBS/INXX2) 

X D: MEAN DIAMETER (IN) 

x L21: DISTANCE FROM CENTERLINE OF DRILL TO STRAIN GAGE (IN) 

X L22: DISTANCE FROM CENTERLINE OF DRILL TO STRAIN GAGE(IN) 



* 


LI i 


DISTANCE FROM DRILL 


TIP TO CYLINDER CENTROID (IN) 


X 


FI : 


X 


AXIS LOAD (LBS) 










* 


F2 i 


Y 


AXIS LOAD (LBS) 








• 


X 


F3 i 


Z 


AXIS LOAD (LBS) 










X 


Ml 1 


MOMENT LOAD ABOUT X 


AXIS 


(IN 


-LBS) 




X 


M2: 


MOMENT LOAD ABOUT Y 


AXIS 


(IN 


-LBS) 




X 


M3: 


MOMENT LOAD ABOUT Z 


AXIS 


(IN 


-LBS) 




X 


EYA21 ; 


NORMAL STRAIN AT 


POINT 


A 


DISTANCE 


L 21 


X 


EYA22 : 


NORMAL STRAIN AT 


POINT 


A 


DISTANCE 


L22 


X 


EYB21 ; 


NORMAL STRAIN AT 


POINT 


B 


DISTANCE 


L21 


X 


EYB22: 


NORMAL STRAIN AT 


POINT 


B 


DISTANCE 


L22 


X 


EYC 


21 : 


NORMAL STRAIN AT 


POINT 


C 


DISTANCE 


L21 


X 


GYZB21 


: SHEAR STRAIN AT 


POINT 


B 


DISTANCE 


L21 


X 


EZB21 : 


NORAML STRAIN AT 


POINT 


3 


DISTANCE 


L21 



x EYZB21: Y-Z PLANE STRAIN ON 45 DEGREE ROSETTE LEG DISTANCE L21 

X X( I ) : CONSTANTS INCLUDE I , Q , RO , D , P I . u , E 

x A(I.J): REPEATED CALCULATIONS IN EQUATIONS 

*XXXXXXXXXXXXX**XXXaXX*XXX*>XX*XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX) 
DOUBLE PRECISION RO , RI , TH . E, POI , PI , G, D, L 21 ,122 , L 1 , XI . X2 , X3 . X4 
+ A1 , A21 , A22 , B21 , B22, EYA21 , EYA22, EYS21 , EYB22, EYC21 , EYZ321 , EZB21 , 
x SYZ321, FI, F2, F3.M1 ,N2,M3, Y1 ,Y2,Y3 
INTEGER J 

XDEFINE KNOWN PARAMETERS 
RO = U . 5 
TH = 0 . 1 
RI=RO-TH 

£=29000000. 0 - . 

P01 = 0 . 29 ■ 

PI=4 . 0*ATAN( 1.0) 

G=E/(2.0X(1 .O+POI)) 

D= ( 2 . 0 *RO ) -TH 
L 21 = 2 . 424 
L22=6 .424 
L1 = 0 . 0 

XCALCULATE CONSTANT TERMS IN EQUATIONS 
X1=PIXTHXDX£ 

X2=4 .OXRO/(PIXEX(ROXX4-RIXX4)) 

X3=8.0x(ROXx3-RIXX3)/(3. 0XGXPIX(THXX2)X(DXX3) ) 

X4=GxPlxTHx(Dxx2) 

XINPUT THE LOADS 
DO 10 J = 1 , 4 

READ(14,x) F1,F2,F3,M1,M2,M3 
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A1=F2/X1 

A21=-(L21XF3+L1XF2+M1)XX2 
A22=-(L22XF3+L1XF2+M1)XX2 
B21=(L21XF1-M3)XX2 
B22=(L22XF1-M3)XX2 
XCALCULATE THE STRAINS 
EYA21 =( A21+A1 ) 

EYA22=( A22+A1 ) 

EYB21=( B21+A1 ) 

EYB22= ( B22+A1 ) 

EYC21=(-A21+A1) 

SY2B21=(-2.0X(L1*F1-M2)/XA)+(X3XF3) 
E2B21=-P0IXEYB21 
EY2B21 = ( EYB21+EZB21+SYZB21 )/2 . 0 
XPRINT LOADS FOR INPUT TO EASYPLOT 



0PEN(UNIT=21 , FILE 3 ’ DATA1 ' , STATUS 3 ’OLD' ) 

WRITEC 21 , 5 ) 

5 FORMAT C T5, » FI • , T13, • F2' , T21 , • F3 ' , T29 , ' Ml ’ , T37 , *M2 • , T45, ' M3 V ) 

WRITEC21 >15) FI , F2 , F3 , Ml ,M2, M3 
XPRINT OUTPUTS OF STRAINS AND LOADS 
15 FORMAT ( 6 FI 0 . 2// ) 

WRITEC 21 / 25 ) 

25 FORMAT CT5, ' EYA21 ’ , T18 , • EYA22 ' ,/) 

WRITEC 21, 35) EYA21,EYA22 
35 FORMAT C 2E12 . A/ ) 

HRITEC21 , 45) 

45 FORMAT CT5, ' EYB21 ' ,T18 , • EYB22* , T30, ' EYZB21 ' , T42, • EZB21 V ) 

WRITEC 21 , 55) EYB21 , EYB22, EYZB21 , EZB21 
55 FORMAT C 4E12 . 4/ ) 

WRITEC 21,65) 

65 F0RMATCT5, 'EYC21'/) 

WRITEC 21 / 75) EYC21 
75 FORMAT C E12 . 4//7 2C ' + *)//) 

XPRINT STRAINS TO FILE FOR CALIBRATION PLT)TS 
EYA21=10U0000XEYA21 
EYA22=1000000X£YA22 
EYB21=1000000XEYB21 
EYB22=1000000XEYB22 
EYZB21=1000000*EYZB21 
EZB21=1000000XEZB21 
EYC21=1000000*EYC21 

WRITEC 20 , X ) EYA21 , EYA22, EYB21 , EYB22, EYZB21 , EZB21 , EYC21 

Y1=EYB22-EYB21 

Y2=EYA21+EYC21 

Y3=EYA22-EYA21 

0PENCUNIT=22,FILE=' DATA2 ' , STATUS = 'OLD' ) 

0PENCUNIT=24,FILE=' DAT A3 STATUS 3 ' OLD ' ) 
0PENCUNIT=25,FILE= , DATA4',STATUS='0LD' ) 

WRITEC 22,85) F1,Y1 
WRITEC 24 , 85 ) F2,Y2 
WRITEC 25, 85)F3,Y3 
85 FORMATC2F10 .2) 

10 CONTINUE 

STOP 
END 



161 



APPENDIX B 



LOAD EQUATIONS FORTRAN COMPUTER CODE 



XXXX> (XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX) 



x LAURENCE, D . A . 

x THIS PROGRAM CALCULATES THE THE LOADS AT THE PENETRATOR DRILL TIP 
X BASED ON MICRO-STRAIN INPUTS. ALSO THE PROGRAM CALCULATES THE LOADS 
x USING THE ASSUMPTIONS THAT F2,M1 , AND M2 ARE ZERO MADE IN THE PREVIOUS 
X WORK. 

X VARIABLES USED ARE: 

X T ; TIME (SEC) 

x TH: RADIAL THICKNESS (IN) 

X RO: OUTSIDE CYLINDER DIAMETER (IN) 

x I: MOMENT OF INERTIA (INXXA) 

X E: YOUNGS MODULUS OF ELASTICITY (LBS/INXX2) 

x POI: POISON RATIO 

X PI 



x 

x 

x 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

< 

X 

X 

X 

X 

X 

X 

X 



Gt SHEAR MODULUS OF ELASTICITY (LBS/INXX2) 

D: MEAN DIAMETER (IN) 

L21 .DISTANCE FROM CENTERLINE OF DRILL TO STRAIN GAGE (IN) 

L22: DISTANCE FROM CENTERLINE OF DRILL TO STRAIN GAGE (IN) 

LI .DISTANCE FROM DRILL TIP TO CYLINDER CENTROID (IN) 

FI: X AXIS LOAD FORCE(LBS) 

F2 : Y AXIS LOAD FORCE(LBS) 

F 3 : Z AXIS LOAD FORCE(LBS) 

Ml: MOMENT LOAD ABOUT X AXIS(IN-LBS) 

M2. MOMENT LOAD ABOUT Y AXIS(IN-LBS) 

M3: MOMENT LOAD ABOUT Z AXIS(IN-LBS) 

EYA21 : NORMAL STRAIN AT POINT A DISTANCE L21 

EYA22: NORMAL STRAIN AT POINT A DISTANCE L22 

EYB21: NORMAL STRAIN AT POINT B DISTANCE L21 

EY322 : NORMAL STRAIN AT POINT B DISTANCE L22 

EYC21 : NORMAL STRAIN AT POHIT C DISTANCE L2I 

SYZE21 : SHEAR STRAIN AT POINT B DISTANCE L21 

EZB21: NORMAL STRAIN AT POINT 3 DISTANCE L21 

EYZB21 : Y-Z PLANE STRAIN ON 45 DEGREE LEG ROSETTE DISTANCE L21 
X( I ) : CONSTANTS INCLUDE I , Q , RO , D, PI , G , E 
A( I , J ) : REPEATED CALCULATIONS IN EQUATIONS 
YF3 : F3 LOAD EQUATION OF PREVIOUS WORK (LBS) 

YM3.M3 LOAD EQUATION OF PREVIOUS WORK (LBS) 

YM2 : M2 LOAD EQUATION OF PREVIOUS WORK (LBS) 



xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx 
DOUBLE PRECISION RO . RI , TH , E, POI , PI , G , D , L21 , L 22 , L 1 , XI , X2 , X3 , 

+ Al, A21, A22, B21, B22.EYA21, EYA22, EYB21 , EYB22, EYC21 , SYZB21 , EZB21 , 

+ EYZB21 , FI , F2, F3 , Ml , M2, M3 , T 
INTEGER J 

XDEFINE KNOWN PARAMETERS 
RO= 0 . 5 



TH = 0 . 1 
RI=RO-TH 
£=29000000. 0 
P01=0 . 29 

PI = 4 . OXATANd .0) 

G=E/(2.0X(1 .O+POI)) 

D=(2.0XR0)-TH 
L21 = 2 . 7 5 
L22 = 6 . 7 5 
Ll=23 . 125 

XCALCULATE CONSTANT TERMS IN EQUATIONS 
Xl=PlxTHxDxE 
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X2=4. OXRO/(PIXE*(RO**4-RI*X4)) 

X3=8. OX!ROxx3-Rlxx3)/!3. 0XGXPIX!THXX2)X(DXX3)) 

X4=GXPIXTHX(DXX2) 

T = 0 . 0 

DO 10 J=1 , 51 

XINPUT MICRO-STRAIN VALUES 

READ( 8 , x ) EYA21 , EYA22, EYB21 , EYB22, EYZB21 , EZB21 , EYC21 
EYA21 = .000001XEYA21 
EYA22= . 000001XEYA22 
EYB21=.000001XEYB21 
EYB22= .000001XEYB22 
EYZB21 = . 000001XEYZB21 
EZB21=.000001XEZB21 
EYC21=.000001XEYC21 
XCALCULATE FORCE AND MOMENT LOADS 

F1=(EYB22-EYB21)/(X2X(L22-L21)) 

F2=X1 X( EYA21 +EYC21 )/2 . 0 

F3 = ! EYA22-EYA21 )/(X2X( -L22+L21 ) ) 

A1=F2/X1 

A21=(L21XF3+L1XF2)XX2 
A22=(L22XF3+L1XF2)XX2 
B21 =-( L21XF1 )XX2 
B22=-( L22XF1 )XX2 
Ml=( (A1-EYA21)/X2)-L21XF3-L1XF2 
M3=( CA1-EYB21)/X2)+L21XF1 
SYZB21 =2 . 0XEYZB21-EZB21-EYB21 
M2=(XAX(SYZB21-X3XF3)/2. 0)+LlXFl 
XPRINT LOADS TO FILE FOR INPUT TO EASYPLOT 
WRITE! 14,2) T,F1,F2,F3,M1,M2,M3 

2 FORMAT C 7F1 0 . 2) 

XCALCULATE LOADS BASED ON PREVIOUS WORK ASSUMPTIONS 
YF3=-EYA21/(X2XL21 ) 

YM3=-EY521/X2 
YM2=X4X(X3XF3-SYZB21 )/2 

XPRINT OLD AND NEW LOADS QUANITIES TO FILE FOR EASYPLOT 
WRITE! 15/ 3) T,F3,YF3,M2,YM2,M3,YM3 

3 FORMAT ! 7 FI 0 . 2 ) 

XPRINTED OUTPUT OF STRAINS AND LOADS 

0PEtl(UNIT = 23,FILE= , DATA5' , STATUS= ' OLD' ) 

WRITE! 23.5)’ 

5 FORMAT ! T5» 'EYA21»,T18, ' EYA22 '/ ) 

WRITE!23/ 15 ) EYA21 , EYA22 
15 . FORMAT ! 2E12 . A/ ) 

WRITE! 23 / 25 ) 

25 FORMAT ! T5, • EYB21 * , T18 , • EYB22 • , T30 , • EYZB21 ' , T42, • EZB21 V ) 

WRITE! 23, 35) EYB21 , EYB22, EYZB21 , EZB21 
35 FORMAT ! 4E1 2 . 4/ ) 

' * WRITE! 23, 45 ) 

45 F0RMAT!T5, • EYC21 ’/) 

WRITE!23,55) EYC21 
55 FORMAT ! E12 . 4/ ) 

WRITE! 23 ,65) 

65 FORMAT! , T , ,T5, , F1',T13, , F2',T21, , F3 , ,T29, 'Ml ' , T37 , 'M2 ' , T45, *M3 V 

+ 48! 'X' )/) 

WRITE! 23, 75) T, FI , F2, F3, Ml , M2 , M3 
75 FORMAT! 7F10.2/72! ' + ' )//) 

T=T+0 .2 

10 CONTINUE 

STOP 
END 
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